From: Akshay S. <ak...@us...> - 2014-02-19 05:11:46
|
This is an automated email from the git hooks/post-receive script. It was generated because a ref change was pushed to the repository containing the project "matlisp". The branch, tensor has been updated via f38e6dde50fbe1552793f8146fa42734d522e9c9 (commit) via 4ae0303bba3df2d7d9b3470181947a0056d72e1b (commit) via 2222db6683c9dbf031cd4db8db5214efe60b6d66 (commit) via 6c30013f4baa53a1b9fba64854c5c1e5cae44809 (commit) via c248fe3323b34374070cb9df9a6d765a85e73b01 (commit) via 17a8a5233aa62740a17e8049835976f7a18e3d26 (commit) via 2e87492c26e3e9f0705efda698f6183d9c1425ea (commit) via 4d63cc7ebed68cf20b1b4e83cbfb6b8815706a4e (commit) via b6f729d172193ff03cf1ba88d1deb1c7634ee11f (commit) via 1c59134bdfcda89a91ce78f8d69836fd3a2628ec (commit) via 7cd35fab7aa468327b733ab1d5037a5e98c55e08 (commit) via e51ecd915cbd2a9222b653d70bda556411616999 (commit) via 983fa49410b5ff5805ef9f63776884fc72015f49 (commit) via 673b1af27a8d2ef318dc02b9b73aa9ce2f758fcc (commit) via ad1dd99286b8c8f0ec1323aaca6911f7f3fd4c99 (commit) via 8a5ade0a47e01bd93e19f72fcfe9691ed00f71cf (commit) via 57618ec426afa04b6553dec19c4c96971c98f5ad (commit) from 7ddfe787e54e485108ff96839495e7a6f0d490c2 (commit) Those revisions listed above that are new to this repository have not appeared on any other notification email; so we list those revisions in full, below. - Log ----------------------------------------------------------------- commit f38e6dde50fbe1552793f8146fa42734d522e9c9 Author: Akshay Srinivasan <aks...@gm...> Date: Tue Feb 18 21:12:52 2014 -0800 Moved all the BLAS things into one folder. Moved gemv, gemm to use define-tensor-method. diff --git a/matlisp.asd b/matlisp.asd index 71c0a25..719f4f8 100644 --- a/matlisp.asd +++ b/matlisp.asd @@ -136,8 +136,8 @@ (:file "symbolic-tensor") (:file "matrix" :depends-on ("numeric")))) - (:module "matlisp-level-1" - :pathname "level-1" + (:module "matlisp-blas" + :pathname "blas" :depends-on ("matlisp-base" "matlisp-classes" "foreign-core") :components ((:file "maker") (:file "copy" @@ -154,31 +154,27 @@ (:file "trans" :depends-on ("scal" "copy")) (:file "sum" - :depends-on ("dot" "copy")))) - (:module "matlisp-level-2" - :pathname "level-2" - :depends-on ("matlisp-base" "matlisp-classes" "foreign-core" "matlisp-level-1") - :components ((:file "gemv"))) - (:module "matlisp-level-3" - :pathname "level-3" - :depends-on ("matlisp-base" "matlisp-classes" "foreign-core" "matlisp-level-1" "matlisp-level-2") - :components ((:file "gemm"))) + :depends-on ("dot" "copy")) + (:file "gemv" + :depends-on ("copy")) + (:file "gemm" + :depends-on ("copy")))) (:module "matlisp-lapack" :pathname "lapack" - :depends-on ("matlisp-base" "matlisp-classes" "matlisp-level-1" "matlisp-level-2" "matlisp-level-3") + :depends-on ("matlisp-base" "matlisp-classes" "matlisp-blas") :components ((:file "lu") (:file "chol") (:file "eig") (:file "least-squares"))) (:module "matlisp-special" :pathname "special" - :depends-on ("matlisp-base" "matlisp-classes" "matlisp-level-1" "matlisp-level-2" "matlisp-level-3") + :depends-on ("matlisp-base" "matlisp-classes" "matlisp-blas") :components ((:file "random") (:file "map") (:file "seq"))) (:module "matlisp-sugar" :pathname "sugar" - :depends-on ("matlisp-base" "matlisp-classes" "matlisp-level-1" "matlisp-level-2" "matlisp-level-3") + :depends-on ("matlisp-base" "matlisp-classes" "matlisp-blas") :components (#+nil (:file "mplusminus") #+nil diff --git a/src/level-1/axpy.lisp b/src/blas/axpy.lisp similarity index 100% rename from src/level-1/axpy.lisp rename to src/blas/axpy.lisp diff --git a/src/level-1/copy.lisp b/src/blas/copy.lisp similarity index 100% rename from src/level-1/copy.lisp rename to src/blas/copy.lisp diff --git a/src/level-1/dot.lisp b/src/blas/dot.lisp similarity index 100% rename from src/level-1/dot.lisp rename to src/blas/dot.lisp diff --git a/src/level-3/gemm.lisp b/src/blas/gemm.lisp similarity index 80% rename from src/level-3/gemm.lisp rename to src/blas/gemm.lisp index aadb836..960b03a 100644 --- a/src/level-3/gemm.lisp +++ b/src/blas/gemm.lisp @@ -117,35 +117,21 @@ (= nc-a nr-b) (= nc-b nc-c)) nil 'tensor-dimension-mismatch)))) -(defmethod gemm! (alpha (A standard-tensor) (B standard-tensor) beta (C standard-tensor) &optional (job :nn)) - (let ((cla (class-name (class-of A))) - (clb (class-name (class-of B))) - (clc (class-name (class-of C)))) - (assert (and (member cla *tensor-type-leaves*) - (member clb *tensor-type-leaves*) - (member clc *tensor-type-leaves*)) - nil 'tensor-abstract-class :tensor-class (list cla clb clc)) - (cond - ((ieql cla clb clc) - (compile-and-eval - `(defmethod gemm! (alpha (A ,cla) (B ,clb) beta (C ,clc) &optional (job :nn)) - (let ((alpha (t/coerce ,(field-type cla) alpha)) - (beta (t/coerce ,(field-type cla) beta))) - (declare (type ,(field-type cla) alpha beta)) - (destructuring-bind (joba jobb) (split-job job) - (declare (type character joba jobb)) - ,(recursive-append - (when (subtypep clc 'blas-numeric-tensor) - `(if (call-fortran? C (t/l3-lb ,clc)) - (with-columnification (,cla ((a joba) (b jobb)) (c)) - (multiple-value-bind (lda opa) (blas-matrix-compatiblep a joba) - (multiple-value-bind (ldb opb) (blas-matrix-compatiblep b jobb) - (t/blas-gemm! ,cla alpha A lda B ldb beta C (or (blas-matrix-compatiblep c #\N) 0) opa opb)))))) - `(t/gemm! ,cla alpha A B beta C joba jobb)))) - C)) - (gemm! alpha A B beta C job)) - (t - (error "Don't know how to apply gemm! to classes ~a." (list cla clb clc)))))) +(define-tensor-method gemm! (alpha (A standard-tensor :input) (B standard-tensor :input) beta (C standard-tensor :output) &optional (job :nn)) + `(let ((alpha (t/coerce ,(field-type (cl a)) alpha)) + (beta (t/coerce ,(field-type (cl a)) beta))) + (declare (type ,(field-type (cl a)) alpha beta)) + (destructuring-bind (joba jobb) (split-job job) + (declare (type character joba jobb)) + ,(recursive-append + (when (subtypep (cl c) 'blas-numeric-tensor) + `(if (call-fortran? C (t/l3-lb ,(cl c))) + (with-columnification (,(cl a) ((a joba) (b jobb)) (c)) + (multiple-value-bind (lda opa) (blas-matrix-compatiblep a joba) + (multiple-value-bind (ldb opb) (blas-matrix-compatiblep b jobb) + (t/blas-gemm! ,(cl a) alpha A lda B ldb beta C (or (blas-matrix-compatiblep c #\N) 0) opa opb)))))) + `(t/gemm! ,(cl a) alpha A B beta C joba jobb)))) + 'C) ;;---------------------------------------------------------------;; (defgeneric gemm (alpha a b beta c &optional job) (:documentation diff --git a/src/level-2/gemv.lisp b/src/blas/gemv.lisp similarity index 75% rename from src/level-2/gemv.lisp rename to src/blas/gemv.lisp index c305481..e468a76 100644 --- a/src/level-2/gemv.lisp +++ b/src/blas/gemv.lisp @@ -93,41 +93,27 @@ (aref (the index-store-vector (dimensions A)) (if (member job '(:t :c)) 1 0)))) nil 'tensor-dimension-mismatch))) -(defmethod gemv! (alpha (A standard-tensor) (x standard-tensor) beta (y standard-tensor) &optional (job :n)) - (let ((clx (class-name (class-of x))) - (cly (class-name (class-of y))) - (cla (class-name (class-of A)))) - (assert (and (member cla *tensor-type-leaves*) - (member clx *tensor-type-leaves*) - (member cly *tensor-type-leaves*)) - nil 'tensor-abstract-class :tensor-class (list cla clx cly)) - (cond - ((ieql clx cly cla) - (compile-and-eval - `(defmethod gemv! (alpha (A ,cla) (x ,clx) beta (y ,cly) &optional (job :n)) - (let ((alpha (t/coerce ,(field-type clx) alpha)) - (beta (t/coerce ,(field-type clx) beta)) - (cjob (aref (symbol-name job) 0))) - (declare (type ,(field-type clx) alpha beta) - (type character cjob)) - ,(recursive-append - (when (subtypep clx 'blas-numeric-tensor) - `(if (call-fortran? A (t/l2-lb ,cla)) - (let ((A-copy (if (blas-matrix-compatiblep A cjob) A - (let ((*default-stride-ordering* :col-major)) - (t/copy! (,cla ,cla) A (t/zeros ,clx (dimensions A))))))) - (multiple-value-bind (lda op maj) (blas-matrix-compatiblep A-copy cjob) - (declare (ignore maj)) - (t/blas-gemv! ,cla alpha A-copy lda - x (aref (the index-store-vector (strides x)) 0) - beta - y (aref (the index-store-vector (strides y)) 0) - op))))) - `(t/gemv! ,cla alpha A x beta y cjob))) - y)) - (gemv! alpha A x beta y job)) - (t - (error "Don't know how to apply gemv! to classes ~a." (list cla clx cly)))))) +(define-tensor-method gemv! (alpha (A standard-tensor :input) (x standard-tensor :input) beta (y standard-tensor :output) &optional (job :n)) + `(let ((alpha (t/coerce ,(field-type (cl x)) alpha)) + (beta (t/coerce ,(field-type (cl x)) beta)) + (cjob (aref (symbol-name job) 0))) + (declare (type ,(field-type (cl x)) alpha beta) + (type character cjob)) + ,(recursive-append + (when (subtypep (cl x) 'blas-numeric-tensor) + `(if (call-fortran? A (t/l2-lb ,(cl a))) + (let ((A-copy (if (blas-matrix-compatiblep A cjob) A + (let ((*default-stride-ordering* :col-major)) + (t/copy! (,(cl a) ,(cl a)) A (t/zeros ,(cl x) (dimensions A))))))) + (multiple-value-bind (lda op maj) (blas-matrix-compatiblep A-copy cjob) + (declare (ignore maj)) + (t/blas-gemv! ,(cl a) alpha A-copy lda + x (aref (the index-store-vector (strides x)) 0) + beta + y (aref (the index-store-vector (strides y)) 0) + op))))) + `(t/gemv! ,(cl a) alpha A x beta y cjob))) + 'y) ;;---------------------------------------------------------------;; (defgeneric gemv (alpha A x beta y &optional job) (:documentation diff --git a/src/level-1/maker.lisp b/src/blas/maker.lisp similarity index 100% rename from src/level-1/maker.lisp rename to src/blas/maker.lisp diff --git a/src/level-1/realimag.lisp b/src/blas/realimag.lisp similarity index 100% rename from src/level-1/realimag.lisp rename to src/blas/realimag.lisp diff --git a/src/level-1/scal.lisp b/src/blas/scal.lisp similarity index 100% rename from src/level-1/scal.lisp rename to src/blas/scal.lisp diff --git a/src/level-1/sum.lisp b/src/blas/sum.lisp similarity index 100% rename from src/level-1/sum.lisp rename to src/blas/sum.lisp diff --git a/src/level-1/swap.lisp b/src/blas/swap.lisp similarity index 100% rename from src/level-1/swap.lisp rename to src/blas/swap.lisp diff --git a/src/level-1/trans.lisp b/src/blas/trans.lisp similarity index 100% rename from src/level-1/trans.lisp rename to src/blas/trans.lisp diff --git a/src/classes/numeric.lisp b/src/classes/numeric.lisp index 018fdc4..2f5c97a 100644 --- a/src/classes/numeric.lisp +++ b/src/classes/numeric.lisp @@ -100,7 +100,7 @@ (imagpart (imagpart element))) (format stream (if (zerop imagpart) "~11,5,,,,,'Eg" - "#C(~0,4,,,,,'Ee, ~0,4,,,,,'Ee)") + "#C(~11,5,,,,,'Eg, ~11,5,,,,,'Eg)") realpart imagpart))) ;; (defleaf complex-tensor (complex-numeric-tensor) ()) diff --git a/src/foreign-core/blas.lisp b/src/foreign-core/blas.lisp index 2e1f57c..202f555 100644 --- a/src/foreign-core/blas.lisp +++ b/src/foreign-core/blas.lisp @@ -28,6 +28,24 @@ (in-package #:matlisp-blas) +;; (defparameter *f77-floats* '(:single-float :double-float :complex-single-float :complex-double-float)) + +;; (defmacro generate-bindings (fname) +;; (let ((defs (parse-fortran-file fname))) +;; `(eval-every +;; ,@(mapcar #'(lambda (x) +;; `(def-fortran-routine ,(first x) ,(second x) +;; ,@(mapcar #'(lambda (y) +;; (let ((type (cadr y)) +;; (var (car y))) +;; (if (and (consp type) (eql (first type) '*) (member (second type) *f77-floats*)) +;; (list var (append type (list :inc (intern (string+ "HEAD-" (symbol-name var)))))) +;; y))) +;; (third x)))) +;; defs)))) + +;; (generate-bindings "/home/neptune/devel/matlisp/blas/blas.f") + (def-fortran-routine daxpy :void " Syntax commit 4ae0303bba3df2d7d9b3470181947a0056d72e1b Author: Akshay Srinivasan <aks...@gm...> Date: Tue Feb 18 20:51:23 2014 -0800 Moved L1 functions to use define-tensor-method. diff --git a/src/base/tensor-template.lisp b/src/base/tensor-template.lisp index 7c22736..ad93e54 100644 --- a/src/base/tensor-template.lisp +++ b/src/base/tensor-template.lisp @@ -147,11 +147,30 @@ (lst (assoc ',(mapcar #'(lambda (x) (if (consp x) (cadr x) t)) args) (cdr (gethash ',name *generated-methods*)) :test #'list-eq))) (assert lst nil "Method table missing from *generated-methods* !") (setf (cdr lst) (list* method (cdr lst)))) - (,name ,@(mapcar #'(lambda (x) (if (consp x) (car x) x)) args))) + (,name ,@(mapcar #'(lambda (x) (if (consp x) (car x) x)) (remove-if #'(lambda (x) (and (not (consp x)) (char= (aref (symbol-name x) 0) #\&))) args)))) ((and (every #'(lambda (,x) (eql ,x (car ,oclasses))) ,oclasses) (or (null ,oclasses) (coerceable? (cclass-max ,iclasses) (car ,oclasses)))) (let* ((clm (or (car ,oclasses) (cclass-max ,iclasses))) ,@(mapcar #'(lambda (x) `(,x (lazy-coerce ,x clm))) inputs)) - (,name ,@(mapcar #'(lambda (x) (if (consp x) (car x) x)) args)))) + (,name ,@(mapcar #'(lambda (x) (if (consp x) (car x) x)) (remove-if #'(lambda (x) (and (not (consp x)) (char= (aref (symbol-name x) 0) #\&))) args))))) (t (error "Don't know how to apply ~a to classes ~a, ~a." ',name ,iclasses ,oclasses))))))))) + + +;; + +;; (defgeneric testg (a)) +;; (define-tensor-method testg ((x standard-tensor :output)) +;; `(t/copy! (t ,(cl x)) 1 x) +;; 'x) + +;; (defgeneric axpy-test (alpha x y)) + +;; (define-tensor-method axpy-test (alpha (x standard-tensor :input) (y standard-tensor :output)) +;; `(let ((alpha (t/coerce ,(field-type (cl x)) alpha))) +;; (declare (type ,(field-type (cl x)) alpha)) +;; ,(recursive-append +;; (when (subtypep (cl x) 'blas-numeric-tensor) +;; `(if-let (strd (and (call-fortran? x (t/l1-lb ,(cl x))) (blas-copyablep x y))) +;; (t/blas-axpy! ,(cl x) alpha x (first strd) y (second strd)))) +;; `(t/axpy! ,(cl x) alpha x y)))) diff --git a/src/level-1/axpy.lisp b/src/level-1/axpy.lisp index ee782f7..2e28522 100644 --- a/src/level-1/axpy.lisp +++ b/src/level-1/axpy.lisp @@ -104,62 +104,25 @@ (assert (lvec-eq (dimensions x) (dimensions y) #'=) nil 'tensor-dimension-mismatch))) -;; - -;; (defgeneric testg (a)) -;; (define-tensor-method testg ((x standard-tensor :output)) -;; `(t/copy! (t ,(cl x)) 1 x) -;; 'x) - -;; (defgeneric axpy-test (alpha x y)) - -;; (define-tensor-method axpy-test (alpha (x standard-tensor :input) (y standard-tensor :output)) -;; `(let ((alpha (t/coerce ,(field-type (cl x)) alpha))) -;; (declare (type ,(field-type (cl x)) alpha)) -;; ,(recursive-append -;; (when (subtypep (cl x) 'blas-numeric-tensor) -;; `(if-let (strd (and (call-fortran? x (t/l1-lb ,(cl x))) (blas-copyablep x y))) -;; (t/blas-axpy! ,(cl x) alpha x (first strd) y (second strd)))) -;; `(t/axpy! ,(cl x) alpha x y)))) - -(defmethod axpy! (alpha (x standard-tensor) (y standard-tensor)) - (let ((clx (class-name (class-of x))) - (cly (class-name (class-of y)))) - (assert (and (member clx *tensor-type-leaves*) - (member cly *tensor-type-leaves*)) - nil 'tensor-abstract-class :tensor-class (list clx cly)) - (cond - ((eq clx cly) - (compile-and-eval - `(defmethod axpy! ((alpha t) (x ,clx) (y ,cly)) - (let ((alpha (t/coerce ,(field-type clx) alpha))) - (declare (type ,(field-type clx) alpha)) - ,(recursive-append - (when (subtypep clx 'blas-numeric-tensor) - `(if-let (strd (and (call-fortran? x (t/l1-lb ,clx)) (blas-copyablep x y))) - (t/blas-axpy! ,clx alpha x (first strd) y (second strd)))) - `(t/axpy! ,clx alpha x y)) - y))) - (axpy! alpha x y)) - (t - (error "Don't know how to apply axpy! to classes ~a, ~a." clx cly))))) - -(defmethod axpy! (alpha (x (eql nil)) (y standard-tensor)) - (let ((cly (class-name (class-of y)))) - (assert (member cly *tensor-type-leaves*) - nil 'tensor-abstract-class :tensor-class cly) - (compile-and-eval - `(defmethod axpy! ((alpha t) (x (eql nil)) (y ,cly)) - (let ((alpha (t/coerce ,(field-type cly) alpha))) - (declare (type ,(field-type cly) alpha)) - ,(recursive-append - (when (subtypep cly 'blas-numeric-tensor) - `(if-let (strd (and (call-fortran? y (t/l1-lb ,cly)) (consecutive-storep y))) - (t/blas-axpy! ,cly alpha nil nil y strd))) - `(t/axpy! ,cly alpha nil y)) - y))) - (axpy! alpha nil y))) - +(define-tensor-method axpy! (alpha (x standard-tensor :input) (y standard-tensor :output)) + `(let ((alpha (t/coerce ,(field-type (cl x)) alpha))) + (declare (type ,(field-type (cl x)) alpha)) + ,(recursive-append + (when (subtypep (cl x) 'blas-numeric-tensor) + `(if-let (strd (and (call-fortran? x (t/l1-lb ,(cl x))) (blas-copyablep x y))) + (t/blas-axpy! ,(cl x) alpha x (first strd) y (second strd)))) + `(t/axpy! ,(cl x) alpha x y)) + y)) + +(define-tensor-method axpy! (alpha (x (eql nil)) (y standard-tensor :output)) + `(let ((alpha (t/coerce ,(field-type (cl y)) alpha))) + (declare (type ,(field-type (cl y)) alpha)) + ,(recursive-append + (when (subtypep (cl y) 'blas-numeric-tensor) + `(if-let (strd (and (call-fortran? y (t/l1-lb ,(cl y))) (consecutive-storep y))) + (t/blas-axpy! ,(cl y) alpha nil nil y strd))) + `(t/axpy! ,(cl y) alpha nil y)) + y)) ;; (defgeneric axpy (alpha x y) (:documentation @@ -186,5 +149,7 @@ (axpy! alpha x (copy y)))) (defmethod axpy (alpha (x standard-tensor) (y (eql nil))) - (let ((tmp (zeros (dimensions x) (class-of x)))) - (axpy! alpha x tmp))) + (axpy! alpha x (zeros (dimensions x) (class-of x)))) + +(defmethod axpy ((alpha complex) (x real-numeric-tensor) (y (eql nil))) + (axpy! alpha x (zeros (dimensions x) 'complex-tensor))) diff --git a/src/level-1/dot.lisp b/src/level-1/dot.lisp index f3d136a..9de26de 100644 --- a/src/level-1/dot.lisp +++ b/src/level-1/dot.lisp @@ -117,46 +117,28 @@ (* (conjugate x) y) (* x y))) -(defmethod dot ((x standard-tensor) (y standard-tensor) &optional (conjugate-p t)) - (let ((clx (class-name (class-of x))) - (cly (class-name (class-of y)))) - (assert (and (member clx *tensor-type-leaves*) - (member cly *tensor-type-leaves*)) - nil 'tensor-abstract-class :tensor-class (list clx cly)) - (cond - ((eq clx cly) - (compile-and-eval - `(defmethod dot ((x ,clx) (y ,cly) &optional (conjugate-p t)) - ,(recursive-append - (when (subtypep clx 'blas-numeric-tensor) - `(if (call-fortran? x (t/l1-lb ,clx)) - (if conjugate-p - (t/blas-dot ,clx x y t) - (t/blas-dot ,clx x y nil)))) - `(if conjugate-p - ;;Please do your checks before coming here. - (t/dot ,clx x y t) - (t/dot ,clx x y nil))))) - (dot x y conjugate-p)) - (t - (error "Don't know how to compute the dot product of ~a , ~a." clx cly))))) +(define-tensor-method dot ((x standard-tensor :input) (y standard-tensor :input) &optional (conjugate-p t)) + (recursive-append + (when (subtypep (cl x) 'blas-numeric-tensor) + `(if (call-fortran? x (t/l1-lb ,(cl x))) + (if conjugate-p + (t/blas-dot ,(cl x) x y t) + (t/blas-dot ,(cl x) x y nil)))) + `(if conjugate-p + ;;Please do your checks before coming here. + (t/dot ,(cl x) x y t) + (t/dot ,(cl x) x y nil)))) -(defmethod dot ((x standard-tensor) (y t) &optional (conjugate-p t)) - (let ((clx (class-name (class-of x)))) - (assert (member clx *tensor-type-leaves*) - nil 'tensor-abstract-class :tensor-class (list clx)) - (compile-and-eval - `(defmethod dot ((x ,clx) (y t) &optional (conjugate-p t)) - (let ((y (t/coerce ,(field-type clx) y))) - (declare (type ,(field-type clx) y)) - ,(recursive-append - (when (subtypep clx 'blas-numeric-tensor) - `(if (call-fortran? x (t/l1-lb ,clx)) - (if conjugate-p - (t/blas-dot ,clx x y t t) - (t/blas-dot ,clx x y nil t)))) - `(if conjugate-p - ;;Please do your checks before coming here. - (t/dot ,clx x y t t) - (t/dot ,clx x y nil t)))))) - (dot x y conjugate-p))) +(define-tensor-method dot ((x standard-tensor :input) (y t) &optional (conjugate-p t)) + `(let ((y (t/coerce ,(field-type (cl x)) y))) + (declare (type ,(field-type (cl x)) y)) + ,(recursive-append + (when (subtypep (cl x) 'blas-numeric-tensor) + `(if (call-fortran? x (t/l1-lb ,(cl x))) + (if conjugate-p + (t/blas-dot ,(cl x) x y t t) + (t/blas-dot ,(cl x) x y nil t)))) + `(if conjugate-p + ;;Please do your checks before coming here. + (t/dot ,(cl x) x y t t) + (t/dot ,(cl x) x y nil t))))) diff --git a/src/level-1/realimag.lisp b/src/level-1/realimag.lisp index 498f93d..c7a7a67 100644 --- a/src/level-1/realimag.lisp +++ b/src/level-1/realimag.lisp @@ -42,14 +42,14 @@ If TENSOR is a scalar, returns its real part. " (etypecase tensor - (real-tensor tensor) - (complex-tensor (let ((*check-after-initializing?* nil)) - (make-instance 'real-tensor - :parent-tensor tensor :store (store tensor) - :dimensions (dimensions tensor) - :strides (map 'index-store-vector #'(lambda (x) (* 2 x)) (the index-store-vector (strides tensor))) - :head (the index-type (* 2 (head tensor)))))) - (number (realpart tensor)))) + ((or real-tensor sreal-tensor) tensor) + ((or complex-tensor scomplex-tensor) (let ((*check-after-initializing?* nil)) + (make-instance (if (typep tensor 'complex-tensor) 'real-tensor 'sreal-tensor) + :parent-tensor tensor :store (store tensor) + :dimensions (dimensions tensor) + :strides (map 'index-store-vector #'(lambda (x) (* 2 x)) (the index-store-vector (strides tensor))) + :head (the index-type (* 2 (head tensor)))))) + (number (cl:realpart tensor)))) (definline tensor-imagpart~ (tensor) " @@ -65,13 +65,13 @@ If TENSOR is a scalar, returns its real part. " (etypecase tensor - (real-tensor tensor) - (complex-tensor (let ((*check-after-initializing?* nil)) - (make-instance 'real-tensor - :parent-tensor tensor :store (store tensor) - :dimensions (dimensions tensor) - :strides (map 'index-store-vector #'(lambda (x) (* 2 x)) (the index-store-vector (strides tensor))) - :head (the index-type (1+ (* 2 (head tensor))))))) + ((or real-tensor sreal-tensor) tensor) + ((or complex-tensor scomplex-tensor) (let ((*check-after-initializing?* nil)) + (make-instance (if (typep tensor 'complex-tensor) 'real-tensor 'sreal-tensor) + :parent-tensor tensor :store (store tensor) + :dimensions (dimensions tensor) + :strides (map 'index-store-vector #'(lambda (x) (* 2 x)) (the index-store-vector (strides tensor))) + :head (the index-type (1+ (* 2 (head tensor))))))) (number (realpart tensor)))) (definline tensor-realpart (tensor) diff --git a/src/level-1/scal.lisp b/src/level-1/scal.lisp index eb89d74..5ef3288 100644 --- a/src/level-1/scal.lisp +++ b/src/level-1/scal.lisp @@ -93,97 +93,24 @@ (assert (very-quickly (lvec-eq (the index-store-vector (dimensions x)) (the index-store-vector (dimensions y)) #'=)) nil 'tensor-dimension-mismatch))) -(defmethod scal! ((x standard-tensor) (y standard-tensor)) - (let ((clx (class-name (class-of x))) - (cly (class-name (class-of y)))) - (assert (and (member clx *tensor-type-leaves*) - (member cly *tensor-type-leaves*)) - nil 'tensor-abstract-class :tensor-class (list clx cly)) - (cond - ((eq clx cly) - (compile-and-eval - `(defmethod scal! ((x ,clx) (y ,cly)) - ,(recursive-append - (when (subtypep clx 'blas-numeric-tensor) - `(if-let (strd (and (call-fortran? x (t/l1-lb ,clx)) (blas-copyablep x y))) - (t/blas-scdi! ,clx x (first strd) y (second strd) t))) - `(t/scdi! ,clx x y :scal? t :numx? nil)) - y)) - (scal! x y)) - (t - (error "Don't know how to apply scal! to classes ~a, ~a." clx cly))))) - -(defmethod scal! ((x t) (y standard-tensor)) - (let ((cly (class-name (class-of y)))) - (assert (member cly *tensor-type-leaves*) - nil 'tensor-abstract-class :tensor-class cly) - (compile-and-eval - `(defmethod scal! ((x t) (y ,cly)) - (let ((x (t/coerce ,(field-type cly) x))) - (declare (type ,(field-type cly) x)) - ,(recursive-append - (when (subtypep cly 'blas-numeric-tensor) - `(if-let (strd (and (call-fortran? y (t/l1-lb ,cly)) (consecutive-storep y))) - (t/blas-scdi! ,cly x nil y strd t))) - `(t/scdi! ,cly x y :scal? t :numx? t)) - y))) - (scal! x y))) - -;;These should've auto-generated. -(defgeneric div! (alpha x) - (:documentation - " - Syntax - ====== - (DIV! alpha x) - - Purpose - ======= - X <- X ./ alpha - - Yes the calling order is twisted. -") - (:method :before ((x standard-tensor) (y standard-tensor)) - (assert (very-quickly (lvec-eq (the index-store-vector (dimensions x)) (the index-store-vector (dimensions y)) #'=)) nil - 'tensor-dimension-mismatch))) - -(defmethod div! ((x standard-tensor) (y standard-tensor)) - (let ((clx (class-name (class-of x))) - (cly (class-name (class-of y)))) - (assert (and (member clx *tensor-type-leaves*) - (member cly *tensor-type-leaves*)) - nil 'tensor-abstract-class :tensor-class (list clx cly)) - (cond - ((eq clx cly) - (compile-and-eval - `(defmethod div! ((x ,clx) (y ,cly)) - ,(recursive-append - (when (subtypep clx 'blas-numeric-tensor) - `(if-let (strd (and (call-fortran? x (t/l1-lb ,clx)) (blas-copyablep x y))) - (t/blas-scdi! ,clx x (first strd) y (second strd) nil))) - `(t/scdi! ,clx x y :scal? nil :numx? nil)) - y)) - (div! x y)) - (t - (error "Don't know how to apply div! to classes ~a, ~a." clx cly))))) - -(defmethod div! ((x t) (y standard-tensor)) - (let ((cly (class-name (class-of y)))) - (assert (member cly *tensor-type-leaves*) - nil 'tensor-abstract-class :tensor-class cly) - (compile-and-eval - `(defmethod div! ((x t) (y ,cly)) - (let ((x (t/coerce ,(field-type cly) x))) - (declare (type ,(field-type cly) x)) - ,(recursive-append - (when (subtypep cly 'blas-numeric-tensor) - `(if-let (strd (and (call-fortran? y (t/l1-lb ,cly)) (consecutive-storep y))) - (t/blas-scdi! ,cly x nil y strd nil))) - `(t/scdi! ,cly x y :scal? nil :numx? t)) - y))) - (div! x y))) +(define-tensor-method scal! ((x standard-tensor :input) (y standard-tensor :output)) + (recursive-append + (when (subtypep (cl x) 'blas-numeric-tensor) + `(if-let (strd (and (call-fortran? x (t/l1-lb ,(cl x))) (blas-copyablep x y))) + (t/blas-scdi! ,(cl x) x (first strd) y (second strd) t))) + `(t/scdi! ,(cl x) x y :scal? t :numx? nil)) + 'y) + +(define-tensor-method scal! ((x t) (y standard-tensor :output)) + `(let ((x (t/coerce ,(field-type (cl y)) x))) + (declare (type ,(field-type (cl y)) x)) + ,(recursive-append + (when (subtypep (cl y) 'blas-numeric-tensor) + `(if-let (strd (and (call-fortran? y (t/l1-lb ,(cl y))) (consecutive-storep y))) + (t/blas-scdi! ,(cl y) x nil y strd t))) + `(t/scdi! ,(cl y) x y :scal? t :numx? t)) + y)) -;; (defgeneric scal (alpha x) (:documentation " @@ -204,7 +131,44 @@ (scal! alpha (copy x))) ;;TODO: There is an issue here when x is not coerceable into the tensor class of alpha (:method ((alpha standard-tensor) (x t)) - (scal! alpha (copy! x (zeros (dimensions alpha) (class-of alpha)))))) + ;;We assume commutation of course. + (scal! x (copy alpha)))) + +;;These should've been auto-generated. +(defgeneric div! (alpha x) + (:documentation + " + Syntax + ====== + (DIV! alpha x) + + Purpose + ======= + X <- X ./ alpha + + Yes the calling order is twisted. +") + (:method :before ((x standard-tensor) (y standard-tensor)) + (assert (very-quickly (lvec-eq (the index-store-vector (dimensions x)) (the index-store-vector (dimensions y)) #'=)) nil + 'tensor-dimension-mismatch))) + +(define-tensor-method div! ((x standard-tensor :input) (y standard-tensor :output)) + (recursive-append + (when (subtypep (cl x) 'blas-numeric-tensor) + `(if-let (strd (and (call-fortran? x (t/l1-lb ,(cl x))) (blas-copyablep x y))) + (t/blas-scdi! ,(cl x) x (first strd) y (second strd) nil))) + `(t/scdi! ,(cl x) x y :scal? nil :numx? nil)) + 'y) + +(define-tensor-method div! ((x t) (y standard-tensor :output)) + `(let ((x (t/coerce ,(field-type (cl y)) x))) + (declare (type ,(field-type (cl y)) x)) + ,(recursive-append + (when (subtypep (cl y) 'blas-numeric-tensor) + `(if-let (strd (and (call-fortran? y (t/l1-lb ,(cl y))) (consecutive-storep y))) + (t/blas-scdi! ,(cl y) x nil y strd nil))) + `(t/scdi! ,(cl y) x y :scal? nil :numx? t)) + y)) (defgeneric div (x y) (:documentation " diff --git a/src/level-1/sum.lisp b/src/level-1/sum.lisp index 7da2868..9e08752 100644 --- a/src/level-1/sum.lisp +++ b/src/level-1/sum.lisp @@ -60,4 +60,3 @@ (declare (ignore axis)) (t/sum ,clx x nil)))) (sum! x y axis))) - diff --git a/src/utilities/string.lisp b/src/utilities/string.lisp index aede8fb..3cb393a 100644 --- a/src/utilities/string.lisp +++ b/src/utilities/string.lisp @@ -30,7 +30,8 @@ returning two values: the string and the number of bytes read." (sb-posix:close fd)) (values data fsize))) - (definline split-seq (test seq &key max-cuts) + (declaim (inline split-seq)) + (defun split-seq (test seq &key max-cuts) "Split a sequence, wherever the given character occurs." (let ((split-list nil) (split-count 0) (deletes nil)) (labels ((left-split (prev i) commit 2222db6683c9dbf031cd4db8db5214efe60b6d66 Author: Akshay Srinivasan <aks...@gm...> Date: Tue Feb 18 19:50:24 2014 -0800 Added more types into the F77 parser. diff --git a/packages.lisp b/packages.lisp index 97b776d..064342e 100644 --- a/packages.lisp +++ b/packages.lisp @@ -122,7 +122,7 @@ (:documentation "Fortran foreign function interface")) (defpackage "MATLISP-BLAS" - (:use #:common-lisp #:matlisp-ffi) + (:use #:common-lisp #:matlisp-ffi #:matlisp-utilities) (:export ;;BLAS Level 1 ;;------------ diff --git a/src/ffi/f77-parser.lisp b/src/ffi/f77-parser.lisp index 0a916be..f26b033 100644 --- a/src/ffi/f77-parser.lisp +++ b/src/ffi/f77-parser.lisp @@ -40,7 +40,7 @@ line continuations. line)) ;; (defparameter *%f77.typemap* - '((("character") :char) + '((("character") :character) (("character*") :string) (("character*1") :string) (("character*6") :string) @@ -52,6 +52,7 @@ line continuations. (("complex*16") :complex-double-float) (("external") (* :void)) (("dimension") nil) + (("logical") :integer) (("none") :void))) ;; (defun %f77.type (line) @@ -98,4 +99,3 @@ line continuations. (when (or (member "function" line :test #'string=) (member "subroutine" line :test #'string=)) (push (parse-procedure line) defns)))))) - commit 6c30013f4baa53a1b9fba64854c5c1e5cae44809 Author: Akshay Srinivasan <aks...@gm...> Date: Mon Feb 17 18:43:12 2014 -0800 Added a fortran definition parser. Fixed a bunch of bugs. diff --git a/packages.lisp b/packages.lisp index d492329..97b776d 100644 --- a/packages.lisp +++ b/packages.lisp @@ -116,7 +116,7 @@ #:foreign-vector #:make-foreign-vector #:foreign-vector-p #:fv-ref #:fv-pointer #:fv-size #:fv-type ;;Interface functions - #:def-fortran-routine + #:def-fortran-routine #:parse-fortran-file #:with-vector-data-addresses ) (:documentation "Fortran foreign function interface")) diff --git a/src/ffi/f77-ffi.lisp b/src/ffi/f77-ffi.lisp index e3a9843..5beee53 100644 --- a/src/ffi/f77-ffi.lisp +++ b/src/ffi/f77-ffi.lisp @@ -150,10 +150,13 @@ ((%f77.callback-type-p type) (let* ((callback-name (second type)) (field-gvar (intern (string+ "*" (symbol-name (gensym (symbol-name var))) "*"))) - (c-callback-code (%f77.def-fortran-callback field-gvar callback-name (third type) (cdddr type)))) + (c-callback-code (%f77.def-fortran-callback field-gvar callback-name (third type) (cdddr type)))) (nconsc callback-code `((defvar ,field-gvar nil) ,@c-callback-code)) (nconsc callback-args `((,field-gvar ,var))) (setq ffi-var `(cffi:callback ,callback-name)))) + ;; + ((and (listp type) (eq (car type) '*) (eq (cadr type) :void)) + (setq ffi-var var)) ;; Can't really enforce "style" when given an array. ;; Complex numbers do not latch onto this case, they ;; are passed by value. @@ -317,7 +320,7 @@ (when (and (%f77.output-p style) (not (eq type :string))) (nconsc return-vars `((,func-var ,ffi-var ,type))))))) - + (let ((retvar (gensym))) `( ,(recursive-append @@ -367,7 +370,7 @@ which copies the vector Y of N double-float's to the vector X. The function name in libblas.a is \"dcopy_\" (by Fortran convention). - (DEF-FORTRAN-ROUTINE DCOPY :void + (DEF-FORTRAN-ROUTINE DCOPY :void (N :integer :input) (X (* :double-float) :output) (INCX :integer :input) @@ -398,13 +401,13 @@ NAME Name of the lisp interface function that will be created. The name of the raw FFI will be derived from NAME via the function MAKE-FFI-NAME. The name of foreign function - (presumable a Fortran Function in an external library) + (presumable a Fortran Function in an external library) will be derived from NAME via MAKE-FORTRAN-NAME. RETURN-TYPE The type of data that will be returned by the external (presumably Fortran) function. - + (MEMBER RETURN-TYPE '(:VOID :INTEGER :SINGLE-FLOAT :DOUBLE-FLOAT :COMPLEX-SINGLE-FLOAT :COMPLEX-DOUBLE-FLOAT)) @@ -453,18 +456,18 @@ as one of the values from the lisp function NAME. - ** Note: In all 3 cases above the input VARIABLE will not be destroyed + ** Note: In all 3 cases above the input VARIABLE will not be destroyed or modified directly, a copy is taken and a pointer of that copy is passed to the (presumably Fortran) external routine. - (OR (* X) :INPUT Array entries are used but not modified. - :STRING) :OUTPUT Array entries need not be initialized on input, + (OR (* X) :INPUT Array entries are used but not modified. + :STRING) :OUTPUT Array entries need not be initialized on input, but will be *modified*. In addition, the array will be returned via the Lisp command VALUES from the lisp function NAME. :INPUT-OUTPUT Like :OUTPUT but initial values on entry may be used. - + The keyword :WORKSPACE is a nickname for :INPUT. The keywords :INPUT-OR-OUTPUT, :WORKSPACE-OUTPUT, :WORKSPACE-OR-OUTPUT are nicknames for :OUTPUT. @@ -546,7 +549,7 @@ (,hidden-var-name ,hack-return-type :output) ,@pars)) (setq hack-return-type :void))) - + `(progn (cffi:defcfun (,fortran-name ,lisp-name) ,(%f77.get-return-type hack-return-type) ,@(%f77.parse-fortran-parameters hack-body)) diff --git a/src/ffi/f77-parser.lisp b/src/ffi/f77-parser.lisp new file mode 100644 index 0000000..0a916be --- /dev/null +++ b/src/ffi/f77-parser.lisp @@ -0,0 +1,101 @@ +(in-package #:matlisp-ffi) + +;;Adapted from cl-blapack. +(defun %f77.tokenize (line) + (declare (type string line)) + (split-seq #'(lambda (c) + (cond + ((member c '(#\Space #\,)) t) + ((member c '(#\( #\))) :keep))) + line)) + +(defun %f77.splitlines (line) + " +Split lines of a Fortran 77 file, whilst removing comments, and taking care of +line continuations. +" + (declare (type string line)) + (split-seq (let ((newline-state 0) + (comment-state nil)) + #'(lambda (c) + (cond + ((member c '(#\Newline)) (setf newline-state 0 + comment-state nil) + :delete) + (newline-state + (incf newline-state) + (cond + ((and (= newline-state 1) (member c '(#\* #\C #\c))) + (setf comment-state t + newline-state nil) + :delete) + ((< newline-state 6) + (if (char= c #\Space) :delete + (progn (setf newline-state nil) + :right))) + ((= newline-state 6) + (progn (setf newline-state nil) + (if (member c '(#\Space #\0)) t :delete))))) + (comment-state :delete)))) + line)) +;; +(defparameter *%f77.typemap* + '((("character") :char) + (("character*") :string) + (("character*1") :string) + (("character*6") :string) + (("integer") :integer) + (("real") :single-float) + (("double" "precision") :double-float) + (("complex") :complex-single-float) + (("double" "complex") :complex-double-float) + (("complex*16") :complex-double-float) + (("external") (* :void)) + (("dimension") nil) + (("none") :void))) +;; +(defun %f77.type (line) + (when-let (type (find line *%f77.typemap* :test #'(lambda (x y) (every #'string= x (car y))))) + (list (cadr type) (nthcdr (length (car type)) line)))) + +(defun parse-fortran-file (fname) + (let ((lines (mapcar #'%f77.tokenize (%f77.splitlines (string-downcase (file->string fname)))))) + (labels ((pointerp (pos line) + (let ((lst (nthcdr (1+ pos) line))) + (when (and (consp lst) (every #'(lambda (x y) (if (eql y t) t (string= x y))) lst '("(" t ")"))) + (cadr lst)))) + (parse-procedure (line) + (let* ((func-name (if (string= (car line) "subroutine") + (cadr line) + (elt line (1+ (position "function" line :test #'string=))))) + (output-type (if (string= (car line) "subroutine") + '("none") + (subseq line 0 (position "function" line :test #'string=)))) + (arguments (mapcar #'(lambda (x) (list x nil nil)) (subseq line (1+ (position "(" line :test #'string=)) (position ")" line :test #'string=))))) + (do ((cline '("") (cond + ((null lines) (error "Cannot find END statement.")) + ((string= (caar lines) "end") (pop lines) nil) + (t (pop lines))))) + ((null cline)) + (when (member cline *%f77.typemap* :test #'(lambda (x y) (every #'string= x (car y)))) + (let ((type (%f77.type cline))) + (if (car type) + (mapcar #'(lambda (x) (when (and (not (second x)) (find (car x) (cadr type) :test #'string=)) + (setf (second x) (car type) + (third x) (pointerp (position (car x) (cadr type) :test #'string=) (cadr type))))) + arguments) + (mapcar #'(lambda (x) (when (and (not (third x)) (find (car x) (cadr type) :test #'string=)) + (setf (third x) (pointerp (position (car x) (cadr type) :test #'string=) (cadr type))))) + arguments))))) + (list (intern (string-upcase func-name)) (car (%f77.type output-type)) (mapcar #'(lambda (x) (list (intern (string-upcase (car x))) + (if (null (third x)) + (second x) + (list '* (second x))))) + arguments))))) + (do ((line '("") (pop lines)) + (defns nil)) + ((null line) defns) + (when (or (member "function" line :test #'string=) + (member "subroutine" line :test #'string=)) + (push (parse-procedure line) defns)))))) + diff --git a/src/utilities/functions.lisp b/src/utilities/functions.lisp index c131fa2..3688c36 100644 --- a/src/utilities/functions.lisp +++ b/src/utilities/functions.lisp @@ -36,7 +36,7 @@ (defun list-eq (a b &optional (test #'eq)) (if (or (atom a) (atom b)) (funcall test a b) - (and (list-eq (car a) (car b)) (list-eq (cdr a) (cdr b) test)))) + (and (list-eq (car a) (car b) test) (list-eq (cdr a) (cdr b) test)))) (defun remmeth (func spls &optional quals) (let ((meth (find-method func quals (mapcar #'(lambda (x) (if (consp x) x (find-class x))) spls) nil))) diff --git a/src/utilities/string.lisp b/src/utilities/string.lisp index 39b6510..aede8fb 100644 --- a/src/utilities/string.lisp +++ b/src/utilities/string.lisp @@ -30,36 +30,52 @@ returning two values: the string and the number of bytes read." (sb-posix:close fd)) (values data fsize))) - (defun split-seq (test seq &key (filter-empty? t) max-cuts from-end?) + (definline split-seq (test seq &key max-cuts) "Split a sequence, wherever the given character occurs." - (if (not from-end?) - (let ((split-list nil) - (split-count 0)) - (loop :for i :from 0 :to (length seq) - :with len := (length seq) - :with prev := 0 - :do (let ((cuts-exceeded? (and max-cuts (>= split-count max-cuts)))) - (when (or (= i len) cuts-exceeded? (funcall test (aref seq i))) - (let* ((str (subseq seq prev (if cuts-exceeded? len i)))) - (when (or cuts-exceeded? (< prev i) (not filter-empty?)) - (incf split-count) - (push str split-list)) - (setf prev (1+ i)))) - (when cuts-exceeded? (return)))) - (values (reverse split-list) (1- split-count))) - (let ((split-list nil) - (split-count 0)) - (loop :for i :from (1- (length seq)) :downto -1 - :with prev := (length seq) - :do (let ((cuts-exceeded? (and max-cuts (>= split-count max-cuts)))) - (when (or (< i 0) cuts-exceeded? (funcall test (aref seq i))) - (let ((str (subseq seq (if cuts-exceeded? 0 (1+ i)) prev))) - (when (or cuts-exceeded? (< (1+ i) prev) (not filter-empty?)) - (incf split-count) - (push str split-list)) - (setf prev i))) - (when cuts-exceeded? (return)))) - (values split-list split-count)))) + (let ((split-list nil) (split-count 0) (deletes nil)) + (labels ((left-split (prev i) + (if (not deletes) + (when (< prev i) + (push (subseq seq prev i) split-list) + (incf split-count)) + (do ((dlst deletes (or (cdr dlst) (cons (1- prev) t))) + (pele i (car dlst)) + (ret nil)) + ((eql dlst t) (progn (setf deletes nil) + (when ret + (push (apply #'string+ ret) split-list) + (incf split-count)))) + (let ((ele (car dlst))) + (when (< (1+ ele) pele) + (push (subseq seq (1+ ele) pele) ret))))))) + (loop :for i :from 0 :to (length seq) + :with len := (length seq) + :with prev := 0 + :do (let ((cmd nil)) + (cond + ((or (= i len) (and max-cuts (>= split-count max-cuts))) + (left-split prev len) + (return)) + ((setf cmd (funcall test (aref seq i))) + (case cmd + (:left + (left-split prev (1+ i)) + (setf prev (1+ i))) + (:right + (left-split prev i) + (setf prev i)) + (:keep + (left-split prev i) + (push (string (aref seq i)) split-list) + (incf split-count) + (setf prev (1+ i))) + (:delete + (push i deletes)) + (t + (left-split prev i) + (setf prev (1+ i))))))))) + (values (nreverse split-list) (1- split-count)))) + ;; (defun splitlines (string) "Split the given string wherever the Carriage-return occurs." commit c248fe3323b34374070cb9df9a6d765a85e73b01 Author: Akshay Srinivasan <aks...@gm...> Date: Sun Feb 16 00:27:37 2014 -0800 Added a new define-tensor-macro, which takes care of coercion and assorted irritations. diff --git a/src/base/tensor-template.lisp b/src/base/tensor-template.lisp index 9543583..7c22736 100644 --- a/src/base/tensor-template.lisp +++ b/src/base/tensor-template.lisp @@ -86,3 +86,72 @@ (assert (null (cdr idx)) nil "given more than one index for linear-store") `(setf (aref (the ,(store-type sym) ,store) (the index-type ,(car idx))) (the ,(field-type sym) ,value))) ;; +;;A helper macro which takes of care of the class checking and stuff. +(defparameter *generated-methods* (make-hash-table)) + +(definline lazy-coerce (x out) + (if (typep x out) x + (copy x out))) + +(defun cclass-max (lst) + (let ((max nil)) + (loop :for ele :in lst + :do (when (or (null max) (and (coerceable? max ele) + (or (not (coerceable? ele max)) + (and (subtypep ele 'blas-numeric-tensor) (subtypep max 'blas-numeric-tensor) + (> (float-digits (coerce 0 (store-element-type ele))) + (float-digits (coerce 0 (store-element-type max)))))))) + (setf max ele))) + max)) + +(defmacro define-tensor-method (name (&rest args) &body body) + (let* ((inputs (mapcar #'car (remove-if-not #'(lambda (x) (and (consp x) (eql (third x) :input))) args))) + (outputs (mapcar #'car (remove-if-not #'(lambda (x) (and (consp x) (eql (third x) :output))) args))) + (iclsym (zipsym inputs)) + (oclsym (zipsym outputs))) + ;; + (multiple-value-bind (val exists?) (gethash name *generated-methods*) + (if exists? + (let ((type-meths (assoc (mapcar #'(lambda (x) (if (consp x) (cadr x) t)) args) (cdr val) :test #'list-eq))) + (if type-meths + (progn + (loop :for ele in (cdr type-meths) + :do (remove-method (symbol-function name) ele)) + (setf (cdr type-meths) nil)) + (setf (cdr val) (list* (list (mapcar #'(lambda (x) (if (consp x) (cadr x) t)) args)) (cdr val))))) + (setf (gethash name *generated-methods*) (list name (list (mapcar #'(lambda (x) (if (consp x) (cadr x) t)) args)))))) + ;; + (with-gensyms (x classes iclasses oclasses) + `(defmethod ,name (,@(mapcar #'(lambda (x) (if (consp x) (subseq x 0 2) x)) args)) + (let* (,@(mapcar #'(lambda (lst) `(,(car lst) (class-name (class-of ,(cadr lst))))) (append iclsym oclsym)) + (,iclasses (list ,@(mapcar #'car iclsym))) + (,oclasses (list ,@(mapcar #'car oclsym))) + (,classes (append ,iclasses ,oclasses))) + (labels ((generate-code (class) + (let ((args (mapcar #'(lambda (x) (if (and (consp x) (member (third x) '(:input :output))) + (list (car x) class) + x)) + '(,@args))) + (ebody (macrolet ((cl (,x) + (let ((slook '(,@(mapcar #'(lambda (x) `(,(cadr x) class)) iclsym) + ,@(mapcar #'(lambda (x) `(,(cadr x) class)) oclsym)))) + (or (cadr (assoc ,x slook)) (error "Can't find class of ~a" ,x))))) + (list ,@body)))) + `(defmethod ,',name (,@args) + ,@ebody)))) + (cond + ((every #'(lambda (,x) (eql ,x (car ,classes))) ,classes) + (assert (member (car ,classes) *tensor-type-leaves*) + nil 'tensor-abstract-class :tensor-class ,classes) + (let* ((method (compile-and-eval (generate-code (car ,classes)))) + (lst (assoc ',(mapcar #'(lambda (x) (if (consp x) (cadr x) t)) args) (cdr (gethash ',name *generated-methods*)) :test #'list-eq))) + (assert lst nil "Method table missing from *generated-methods* !") + (setf (cdr lst) (list* method (cdr lst)))) + (,name ,@(mapcar #'(lambda (x) (if (consp x) (car x) x)) args))) + ((and (every #'(lambda (,x) (eql ,x (car ,oclasses))) ,oclasses) + (or (null ,oclasses) (coerceable? (cclass-max ,iclasses) (car ,oclasses)))) + (let* ((clm (or (car ,oclasses) (cclass-max ,iclasses))) + ,@(mapcar #'(lambda (x) `(,x (lazy-coerce ,x clm))) inputs)) + (,name ,@(mapcar #'(lambda (x) (if (consp x) (car x) x)) args)))) + (t + (error "Don't know how to apply ~a to classes ~a, ~a." ',name ,iclasses ,oclasses))))))))) diff --git a/src/foreign-core/blas.lisp b/src/foreign-core/blas.lisp index d7dd57b..2e1f57c 100644 --- a/src/foreign-core/blas.lisp +++ b/src/foreign-core/blas.lisp @@ -28,7 +28,6 @@ (in-package #:matlisp-blas) - (def-fortran-routine daxpy :void " Syntax @@ -72,6 +71,16 @@ (incy :integer :input) ) +(def-fortran-routine saxpy :void + (n :integer :input) + (da :single-float :input) + (dx (* :single-float :inc head-x)) + (incx :integer :input) + (dy (* :single-float :inc head-y) :output) + (incy :integer :input) +) +;; + (def-fortran-routine dcopy :void " Syntax @@ -113,6 +122,15 @@ (incy :integer :input) ) +(def-fortran-routine scopy :void + (n :integer :input) + (dx (* :single-float :inc head-x)) + (incx :integer :input) + (dy (* :single-float :inc head-y) :output) + (incy :integer :input) +) +;; + (def-fortran-routine drot :void " Applies a plane rotation. @@ -245,6 +263,16 @@ (incy :integer :input) ) +(def-fortran-routine caxpy :void + (n :integer :input) + (za :complex-single-float) + (zx (* :complex-single-float :inc head-x)) + (incx :integer :input) + (zy (* :complex-single-float :inc head-y) :output) + (incy :integer :input) +) + +;; (def-fortran-routine zcopy :void " Syntax @@ -287,6 +315,15 @@ (incy :integer :input) ) +(def-fortran-routine ccopy :void + (n :integer :input) + (zx (* :complex-single-float :inc head-x)) + (incx :integer :input) + (zy (* :complex-single-float :inc head-y) :output) + (incy :integer :input) +) + +;; (def-fortran-routine zdscal :void " Syntax diff --git a/src/level-1/axpy.lisp b/src/level-1/axpy.lisp index 6e6f155..ee782f7 100644 --- a/src/level-1/axpy.lisp +++ b/src/level-1/axpy.lisp @@ -30,9 +30,13 @@ (deft/generic (t/blas-axpy-func #'subfieldp) sym ()) (deft/method t/blas-axpy-func (sym real-tensor) () 'daxpy) +(deft/method t/blas-axpy-func (sym sreal-tensor) () + 'saxpy) (deft/method t/blas-axpy-func (sym complex-tensor) () 'zaxpy) -;; +(deft/method t/blas-axpy-func (sym scomplex-tensor) () + 'caxpy) +;; (deft/generic (t/blas-axpy! #'subtypep) sym (a x st-x y st-y)) (deft/method t/blas-axpy! (sym blas-numeric-tensor) (a x st-x y st-y) (let ((apy? (null x))) @@ -100,6 +104,24 @@ (assert (lvec-eq (dimensions x) (dimensions y) #'=) nil 'tensor-dimension-mismatch))) +;; + +;; (defgeneric testg (a)) +;; (define-tensor-method testg ((x standard-tensor :output)) +;; `(t/copy! (t ,(cl x)) 1 x) +;; 'x) + +;; (defgeneric axpy-test (alpha x y)) + +;; (define-tensor-method axpy-test (alpha (x standard-tensor :input) (y standard-tensor :output)) +;; `(let ((alpha (t/coerce ,(field-type (cl x)) alpha))) +;; (declare (type ,(field-type (cl x)) alpha)) +;; ,(recursive-append +;; (when (subtypep (cl x) 'blas-numeric-tensor) +;; `(if-let (strd (and (call-fortran? x (t/l1-lb ,(cl x))) (blas-copyablep x y))) +;; (t/blas-axpy! ,(cl x) alpha x (first strd) y (second strd)))) +;; `(t/axpy! ,(cl x) alpha x y)))) + (defmethod axpy! (alpha (x standard-tensor) (y standard-tensor)) (let ((clx (class-name (class-of x))) (cly (class-name (class-of y)))) diff --git a/src/level-1/copy.lisp b/src/level-1/copy.lisp index 736a40c..12c6b41 100644 --- a/src/level-1/copy.lisp +++ b/src/level-1/copy.lisp @@ -30,8 +30,12 @@ (deft/generic (t/blas-copy-func #'subfieldp) sym ()) (deft/method t/blas-copy-func (sym real-tensor) () 'dcopy) +(deft/method t/blas-copy-func (sym sreal-tensor) () + 'scopy) (deft/method t/blas-copy-func (sym complex-tensor) () 'zcopy) +(deft/method t/blas-copy-func (sym scomplex-tensor) () + 'ccopy) ;; (deft/generic (t/blas-copy! #'subtypep) sym (x st-x y st-y)) (deft/method t/blas-copy! (sym blas-numeric-tensor) (x st-x y st-y) commit 17a8a5233aa62740a17e8049835976f7a18e3d26 Author: Akshay Srinivasan <aks...@gm...> Date: Sun Feb 16 00:25:48 2014 -0800 Updated packages.lisp. diff --git a/packages.lisp b/packages.lisp index 21c0a96..d492329 100644 --- a/packages.lisp +++ b/packages.lisp @@ -77,7 +77,7 @@ #:cut-cons-chain! #:slot-values #:remmeth #:recursive-append #:unquote-args #:flatten - #:format-to-string #:string+ + #:format-to-string #:string+ #:file->string #:split-seq #:splitlines #:linear-array-type #:list-dimensions #:lvec-foldl #:lvec-foldr #:lvec-max #:lvec-min #:lvec-eq @@ -129,9 +129,13 @@ ;;Real-double #:ddot #:dnrm2 #:dasum #:dscal #:daxpy #:drot #:dswap #:dcopy #:idamax + ;;Real-single + #:saxpy #:scopy ;;Complex-double #:zdotc #:zdotu #:zdscal #:zscal #:zswap #:zcopy #:zaxpy #:dcabs1 #:dzasum #:dznrm2 #:izamax + ;;Complex-single + #:caxpy #:ccopy ;;BLAS Level 2 ;;------------ ;;Real-double commit 2e87492c26e3e9f0705efda698f6183d9c1425ea Author: Akshay Srinivasan <aks...@gm...> Date: Sun Feb 16 00:24:41 2014 -0800 Moved string functions into utilities. diff --git a/src/reader/loadsave.lisp b/src/reader/loadsave.lisp index c29cc64..dc9c986 100644 --- a/src/reader/loadsave.lisp +++ b/src/reader/loadsave.lisp @@ -1,45 +1,5 @@ (in-package #:matlisp) -#+(not :sbcl) -(defun file->string (path) - "Sucks up an entire file from PATH into a freshly-allocated string, -returning two values: the string and the number of bytes read." - (declare (optimize (safety 0) (speed 3))) - (with-open-file (s path :external-format :iso8859-1) - (let* ((len (file-length s)) - (data (make-array len :element-type 'standard-char))) - (values data (read-sequence data s))))) - -#+sbcl -(defun file->string (path) -"Sucks up an entire file from PATH into a freshly-allocated string, -returning two values: the string and the number of bytes read." - (let* ((fsize (with-open-file (s path) - (file-length s))) - (data (make-array fsize :element-type 'standard-char)) - (fd (sb-posix:open path 0))) - (unwind-protect (sb-posix:read fd (sb-sys:vector-sap data) fsize) - (sb-posix:close fd)) - (values data fsize))) -;; -(definline split-seq (test seq &optional (filter-empty? t)) - "Split a string, wherever the given character occurs." - (loop :for i :from (1- (length seq)) :downto -1 - :with prev := (length seq) - :with split-list := nil - :with split-count := 0 - :do (when (or (< i 0) (funcall test (aref seq i))) - (let ((str (subseq seq (1+ i) prev))) - (when (or (< (1+ i) prev) (not filter-empty?)) - (incf split-count) - (push str split-list)) - (setf prev i))) - :finally (return (values split-list split-count)))) -;; -(defun splitlines (string) - "Split the given string wherever the Carriage-return occurs." - (split-seq #'(lambda (x) (or (char= x #\Newline) (char= x #\Return))) string)) - ;; ;; (defmacro apply* ((&rest funcl) expr) ;; (let ((syms (zip (mapcar #'gensym funcl) funcl))) diff --git a/src/utilities/string.lisp b/src/utilities/string.lisp index 1833467..39b6510 100644 --- a/src/utilities/string.lisp +++ b/src/utilities/string.lisp @@ -7,4 +7,62 @@ (defun format-to-string (fmt &rest args) (apply #'format (append (list nil fmt) args))) + + #+(not :sbcl) + (defun file->string (path) + "Sucks up an entire file from PATH into a freshly-allocated string, +returning two values: the string and the number of bytes read." + (declare (optimize (safety 0) (speed 3))) + (with-open-file (s path :external-format :iso8859-1) + (let* ((len (file-length s)) + (data (make-array len :element-type 'standard-char))) + (values data (read-sequence data s))))) + + #+sbcl + (defun file->string (path) + "Sucks up an entire file from PATH into a freshly-allocated string, +returning two values: the string and the number of bytes read." + (let* ((fsize (with-open-file (s path) + (file-length s))) + (data (make-array fsize :element-type 'standard-char)) + (fd (sb-posix:open path 0))) + (unwind-protect (sb-posix:read fd (sb-sys:vector-sap data) fsize) + (sb-posix:close fd)) + (values data fsize))) + + (defun split-seq (test seq &key (filter-empty? t) max-cuts from-end?) + "Split a sequence, wherever the given character occurs." + (if (not from-end?) + (let ((split-list nil) + (split-count 0)) + (loop :for i :from 0 :to (length seq) + :with len := (length seq) + :with prev := 0 + :do (let ((cuts-exceeded? (and max-cuts (>= split-count max-cuts)))) + (when (or (= i len) cuts-exceeded? (funcall test (aref seq i))) + (let* ((str (subseq seq prev (if cuts-exceeded? len i)))) + (when (or cuts-exceeded? (< prev i) (not filter-empty?)) + (incf split-count) + (push str split-list)) + (setf prev (1+ i)))) + (when cuts-exceeded? (return)))) + (values (reverse split-list) (1- split-count))) + (let ((split-list nil) + (split-count 0)) + (loop :for i :from (1- (length seq)) :downto -1 + :with prev := (length seq) + :do (let ((cuts-exceeded? (and max-cuts (>= split-count max-cuts)))) + (when (or (< i 0) cuts-exceeded? (funcall test (aref seq i))) + (let ((str (subseq seq (if cuts-exceeded? 0 (1+ i)) prev))) + (when (or cuts-exceeded? (< (1+ i) prev) (not filter-empty?)) + (incf split-count) + (push str split-list)) + (setf prev i))) + (when cuts-exceeded? (return)))) + (values split-list split-count)))) + ;; + (defun splitlines (string) + "Split the given string wherever the Carriage-return occurs." + (split-seq #'(lambda (x) (or (char= x #\Newline) (char= x #\Return))) string)) + ) commit 4d63cc7ebed68cf20b1b4e83cbfb6b8815706a4e Author: Akshay Srinivasan <aks...@gm...> Date: Wed Feb 5 04:50:02 2014 -0800 Saving changes to sparse-tensor.lisp diff --git a/src/base/sparse-tensor.lisp b/src/base/sparse-tensor.lisp index a46c843..4119a3f 100644 --- a/src/base/sparse-tensor.lisp +++ b/src/base/sparse-tensor.lisp @@ -1,6 +1,9 @@ (in-package :matlisp) -;; +;;One may to do better than a Hash-table for this. +(defparameter *default-sparsity* 1/1000) +(defparameter *max-size* 10000) + (defclass coordinate-sparse-tensor (sparse-tensor) ((strides :initarg :strides :reader strides :type index-store-vector :documentation "Strides for accesing elements of the tensor."))) @@ -25,17 +28,17 @@ (unless (t/f= ,(field-type sym) ,val (t/fid+ ,(field-type sym))) (setf (gethash ,(car idx) ,store) (the ,(field-type sym) ,value)))))) +(deft/method t/store-type (sym coordinate-sparse-tensor) (&optional (size '*)) + 'hash-table) + (deft/method t/store-size (sym coordinate-sparse-tensor) (ele) `(hash-table-count ,ele)) (deft/method t/store-type (sym coordinate-sparse-tensor) (&optional (size '*)) 'hash-table) -(defparameter *default-sparsity* 1/1000) -(defparameter *max-size* 10000) - (deft/method t/compute-store-size (sym coordinate-sparse-tensor) (size) - `(max (min sb-impl::+min-hash-table-size+ (ceiling (/ ,size *default-sparsity*))) *max-size*)) + `(max (min sb-impl::+min-hash-table-size+ (ceiling (/ ,size *default-sparsity*))) *max-sparse-size*)) (defmethod head ((tensor coordinate-sparse-tensor)) 0) @@ -44,7 +47,6 @@ (deft/method t/field-type (sym real-sparse-tensor) () 'double-float) - ;; (defmethod ref ((tensor coordinate-sparse-tensor) &rest subscripts) (let ((clname (class-name (class-of tensor)))) @@ -65,61 +67,47 @@ (sto (store tensor))) (t/store-set ,clname (t/coerce ,(field-type clname) value) sto idx) (t/store-ref ,clname sto idx)))) - (setf (ref tensor subscripts) value))) + (setf (ref tensor (if (numberp (car subscripts)) subscripts (car subscripts))) value))) ;; (defclass compressed-sparse-matrix (sparse-tensor) - ((indices :initarg :strides :reader indices :type index-store-vector - :documentation "Strides for accesing elements of the tensor.") - (index-position :initarg :strides :reader index-position :type index-store-vector - :documentation "Strides for accesing elements of the tensor."))) - -(deft/method t/store-allocator (sym compressed-sparse-matrix) (size &optional initial-element) - (with-gensyms (sitm size-sym arr idx init) - (let ((type (macroexpand-1 `(t/store-element-type ,sym)))) - `(let*-typed ((,size-sym (t/compute-store-size ,sym (let ((,sitm ,size)) - (etypecase ,sitm - (index-type ,sitm) - (index-store-vector (lvec-foldr #'* (the index-store-vector ,sitm))) - (cons (reduce #'* ,sitm)))))) - ,@(when initial-element `((,init ,initial-element :type ,(field-type sym)))) - (,arr (make-array ,size-sym :element-type ',type :initial-element ,(if (subtypep type 'number) `(t/f... [truncated message content] |