[lisplab-cvs] r84 - doc/manual src/fft src/linalg src/matlisp
Jørn Inge Vestgården
jivestgarden at common-lisp.net
Sun Aug 16 09:46:41 UTC 2009
Author: jivestgarden
Date: Sun Aug 16 05:46:41 2009
New Revision: 84
Log:
cleaning and error handling
Modified:
Makefile
README
TODO
doc/manual/Makefile
lisplab.asd
package.lisp
src/fft/level3-fft-fftw.lisp
src/fft/level3-fft-zge.lisp
src/linalg/level3-linalg-interface.lisp
src/matlisp/geev.lisp
src/matlisp/inv.lisp
src/matlisp/mul.lisp
start.lisp
Modified: Makefile
==============================================================================
--- Makefile (original)
+++ Makefile Sun Aug 16 05:46:41 2009
@@ -1,7 +1,12 @@
+# Makefile for admin tasks
+.PHONY: first, manual, touch, lispclean, clean, distclean
first:
echo "Please specify target."
+manual:
+ make -C"doc/manual" all
+
touch:
touch system/lisplab.asd
Modified: README
==============================================================================
--- README (original)
+++ README Sun Aug 16 05:46:41 2009
@@ -1,56 +1,45 @@
-INTRODUCTION
-Lisplab is a mathematics library for common lisp, released under GPL.
+============================================================
+LISPLAB - a mathematics library for Common Lisp
+============================================================
+Lisplab is a mathematics library for common lisp, released
+under the GNU General Public License (GPL), except files
+that state something else.
-Directory structure:
- src/ Lisplab source code
- src/matliap From the Matlisp project. Rewritten for Lisplab matrices.
- shared/ Unmodified code from other projects.
- doc/ Documentation.
+Homepage: http://common-lisp.net/projects/lisplab/
+Manual: doc/manual/lisplab.texi
-STATUS
-The project is not in finished state and names will change without warning.
-However, the basic matrix code is a level where it is useful for doing
-mathematical modelling.
-
-NAMING CONVENTION
-The files including "interface" in the name defines generic functions only.
-The files including "generic" in the name defines only unspecialized methods.
-The levels in the names of source code is dependency levels:
-level 0:
- Methods that work on non-index objects,
- but the methods here are also typically reimplemented
- for matrices.
- Exampels: copy, .+, .*, ...
-level 1:
- Basic matrix/tensor methods for indexing, element reference
- and dimensionality. In order to make a new kind of matrix, all
- mathods on this level must be reimplmented.
- Examples: mref, dim, rows, cols, new, ...
-level 2:
- Basic functionality related to matrices.
- Examples: mmax, mmap, diag, ...
-level 3:
- Linear algebra and anything else based on the matrices.
- Examples: minv, m*, ...
-
+============================================================
INSTALLING
-Lisplab is asdf-installable. It has only been tested on SBCL on Linux,
-but should be fairly portable to other platforms, as soon as some minor
-dependencies of the package sb-ext are resolved.
-
-The Matlisp linear algebra depends on externally libraries and these must
-be specified in the variable asdf:*lisplab-external-libraries* before loading,
-as seen in start.lisp. The order of the libraries matter and Blas must be
-before Lapack.
-
-On Linux the Blas and Lapack libraries can often be installed by the operating
-system package system. In Ubuntu, to get the Atlas build of Blas/Lapack, type
- % aptitude install libatlas3gf-base
-
-
-GOOD TO KNOW
-Lisplab only works with double floats and single-floats should not be used.
-To ensure this, use
+============================================================
+Lisplab is asdf-installable. It has only been tested on SBCL
+on Linux, but should be fairly portable to other platforms.
+
+Lisplab uses BLAS, LAPACK, and FFTW, and the path to these
+libraries must be set before loading. See file "start.lisp"
+for how.
+
+Lisplab has only been run with double-floats, so
+your ".sbclrc" should include
(setf *READ-DEFAULT-FLOAT-FORMAT* 'double-float)
+
+
+============================================================
+Organization of the source code
+============================================================
+ src/ Lisplab source code.
+ shared/ Code from Slatec and Quadpack.
+ doc/ Documentation.
+
+
+============================================================
+STATUS
+============================================================
+Think of it as late alpha cycle. The Library is rather extensive,
+the API is fairly stable, and the implementation has OK structure.
+It also has fairly good documentation and a manual.
+But it is poorly tested.
+
+(In general, Lisplab needs more users and more developers.)
+
Modified: TODO
==============================================================================
--- TODO (original)
+++ TODO Sun Aug 16 05:46:41 2009
@@ -1,21 +1,10 @@
TODO
-o About arrays. Handle them by making wrappers in the matrix
- hierarcy. This will be rather slow but saves work.
- And in this way it is possible to avoid the non-speialized methods.
-o Matrices with general element type.
-o Some way to just extend the exact input matrix type
-o Documentation.
o Check if there is some non-threadsafe code in the fortran interface,
i.e. array pre-alocation for the workspaces.
-o Find out how to dynamically switch between common-lisp specialized
- blas-arrays and fortan specialized blas-arrays. Currently this
- is a mess.
-o Test code.
+o Testing.
o Error handling.
-o Added spcialized matrix types, in an ordered way.
-o Package structure.
+o More matrix types. (Sparse matrices)
Extensions:
-o Symbolic maniputions, similar to Ginac in C++.
-o Threaded and paralell execution. Use CUDA?
-
\ No newline at end of file
+o Symbolic manipulations - similar to C++s GINAC.
+o Threaded and paralell execution. Use CUDA, SMP, or just threads?.
Modified: doc/manual/Makefile
==============================================================================
--- doc/manual/Makefile (original)
+++ doc/manual/Makefile Sun Aug 16 05:46:41 2009
@@ -4,10 +4,10 @@
all: html pdf
-html: $(source)
+html:
makeinfo --html $(source)
-pdf: $(source)
+pdf:
texi2pdf $(source)
clean:
Modified: lisplab.asd
==============================================================================
--- lisplab.asd (original)
+++ lisplab.asd Sun Aug 16 05:46:41 2009
@@ -2,23 +2,22 @@
(in-package :cl-user)
-(defvar *lisplab-libblas-path* nil "Path to blas shared object file")
-(defvar *lisplab-liblapack-path* nil "Path to lapack shared object file")
-(defvar *lisplab-fftw-path* nil "Path to fftw shared object file")
+(defvar *lisplab-libblas-path* nil "Path to BLAS shared object file.")
+(defvar *lisplab-liblapack-path* nil "Path to LAPACK shared object file.")
+(defvar *lisplab-libfftw-path* nil "Path to FFTW shared object file.")
(defpackage :asdf-lisplab (:use :asdf :cl))
(in-package :asdf-lisplab)
(defun load-lisplab-lib (name)
- (sb-alien:load-shared-object name))
+ (when name
+ (sb-alien:load-shared-object name)))
+
+(defun explain-lisplab-lib (name path)
+ (format t "Loads ~A. Path ~a" name path))
-(defsystem :lisplab
- ;; Default system, without external libs
- :depends-on
- (:lisplab-base
- :quadpack))
-(defsystem :lisplab-all
+(defsystem :lisplab
;; Default system, without all libs
:depends-on
(:lisplab-base
@@ -136,15 +135,17 @@
(load-lisplab-lib
cl-user::*lisplab-libblas-path*))
:explain (asdf:load-op :after (op c)
- (format t "Loads alien object <~A>"
- cl-user::*lisplab-libblas-path*)))
+ (explain-lisplab-lib
+ "BLAS"
+ cl-user::*lisplab-libblas-path*)))
(:module :lapack-libs
:perform (asdf:load-op :after (op c)
(load-lisplab-lib
cl-user::*lisplab-liblapack-path*))
:explain (asdf:load-op :after (op c)
- (format t "Loads alien object <~A>"
- cl-user::*lisplab-liblapack-path*)))
+ (explain-lisplab-lib
+ "LAPACK"
+ cl-user::*lisplab-liblapack-path*)))
(:file "f77-package")
(:file "ffi-sbcl")
(:file "blas")
@@ -171,8 +172,9 @@
(load-lisplab-lib
cl-user::*lisplab-libfftw-path*))
:explain (asdf:load-op :after (op c)
- (format t "Loads alien object <~A>"
- cl-user::*lisplab-libfftw-path*)))
+ (explain-lisplab-lib
+ "FFTW"
+ cl-user::*lisplab-libfftw-path*)))
(:file "fftw-ffi")
(:file "level3-fft-fftw")))))
Modified: package.lisp
==============================================================================
--- package.lisp (original)
+++ package.lisp Sun Aug 16 05:46:41 2009
@@ -17,7 +17,13 @@
;;; with this program; if not, write to the Free Software Foundation, Inc.,
;;; 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
-;;; TODO split this out
+
+;;; It could be possible to split the big Lisplab package
+;;; in smaller packages as indicated by the comments, but
+;;; this would be more work and I see litle gain from it.
+;;; As it is now, the asdf system tells the dependencies, not
+;;; the pacakges. Think its OK.
+
(defpackage "FORTRAN-FFI-ACCESSORS"
(:use "COMMON-LISP" "SB-ALIEN" "SB-C")
@@ -29,7 +35,17 @@
(defpackage "LISPLAB"
(:use "COMMON-LISP")
- (:export
+ (:nicknames "LL")
+ (:documentation
+ "Lisplab is mathematics library released under the
+GNU General Public License (GPL).
+
+Lisplab contains mathematical functions, matrices,
+linear algebra, Fast Fourier Transform,
+diff-solvers, and a lot more.
+
+Lisplab provides high level interfaces to BLAS, LAPACK and FFTW.")
+ (:export
;; Numerical constants
"%I"
"%E"
@@ -195,15 +211,3 @@
"FFT-SHIFT"
"IFFT-SHIFT"
))
-
-(defpackage "LISPLAB-BLAS"
- (:use "COMMON-LISP" "LISPLAB")
- (:documentation "Mathematics and linear algebra library."))
-
-
-
-
-
-
-
-
Modified: src/fft/level3-fft-fftw.lisp
==============================================================================
--- src/fft/level3-fft-fftw.lisp (original)
+++ src/fft/level3-fft-fftw.lisp Sun Aug 16 05:46:41 2009
@@ -16,13 +16,16 @@
;;; with this program; if not, write to the Free Software Foundation, Inc.,
;;; 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+;;; There is a sligh missuse of notation here, since the method
+;;; specialize on matrix-blas-zge, rather than matrix-fftw-zge, which
+;;; would have a more correct name.
+
(in-package :lisplab)
(defun fft1!-forward-or-backward (x direction)
(let* ((rows (rows x))
(cols (cols x))
(store (matrix-store x)))
- (print "allois")
(dotimes (i cols)
(fftw-ffi:fftw-fft1
rows
@@ -35,27 +38,37 @@
x)
(defmethod fft1! ((x matrix-blas-zge))
- (fft1!-forward-or-backward x fftw-ffi:+fftw-forward+))
+ (if cl-user::*lisplab-libfftw-path*
+ (fft1!-forward-or-backward x fftw-ffi:+fftw-forward+)
+ (call-next-method)))
(defmethod ifft1! ((x matrix-blas-zge))
- (fft1!-forward-or-backward x fftw-ffi:+fftw-backward+))
+ (if cl-user::*lisplab-libfftw-path*
+ (fft1!-forward-or-backward x fftw-ffi:+fftw-backward+)
+ (call-next-method)))
(defmethod fft2! ((x matrix-blas-zge))
- (fftw-ffi:fftw-fft2
- (rows x)
- (cols x)
- (matrix-store x)
- (matrix-store x)
- fftw-ffi:+fftw-forward+
- fftw-ffi:+FFTW-ESTIMATE+)
- x)
+ (if cl-user::*lisplab-libfftw-path*
+ (progn
+ (fftw-ffi:fftw-fft2
+ (rows x)
+ (cols x)
+ (matrix-store x)
+ (matrix-store x)
+ fftw-ffi:+fftw-forward+
+ fftw-ffi:+FFTW-ESTIMATE+)
+ x)
+ (call-next-method)))
(defmethod ifft2! ((x matrix-blas-zge))
- (fftw-ffi:fftw-fft2
- (rows x)
- (cols x)
- (matrix-store x)
- (matrix-store x)
- fftw-ffi:+fftw-backward+
- fftw-ffi:+FFTW-ESTIMATE+)
- x)
+ (if cl-user::*lisplab-libfftw-path*
+ (progn
+ (fftw-ffi:fftw-fft2
+ (rows x)
+ (cols x)
+ (matrix-store x)
+ (matrix-store x)
+ fftw-ffi:+fftw-backward+
+ fftw-ffi:+FFTW-ESTIMATE+)
+ x)
+ (call-next-method)))
Modified: src/fft/level3-fft-zge.lisp
==============================================================================
--- src/fft/level3-fft-zge.lisp (original)
+++ src/fft/level3-fft-zge.lisp Sun Aug 16 05:46:41 2009
@@ -21,25 +21,37 @@
(in-package :lisplab)
+
+
+
;;;; The implementing methods
+
+
(defmethod fft1! ((x matrix-lisp-zge))
+ (assert (= 1 (logcount (rows x))))
(dotimes (i (cols x))
(fft-radix-2-blas-complex-store! :f (matrix-store x) (rows x) (* (rows x) i) 1))
x)
(defmethod ifft1! ((x matrix-lisp-zge))
+ (assert (= 1 (logcount (rows x))))
(dotimes (i (cols x))
(fft-radix-2-blas-complex-store! :r (matrix-store x) (rows x) (* (rows x) i) 1))
x)
(defmethod fft2! ((x matrix-lisp-zge))
+ (assert (and (= 1 (logcount (rows x)))
+ (= 1 (logcount (cols x)))))
(fft1! x)
(dotimes (i (rows x))
(fft-radix-2-blas-complex-store! :f (matrix-store x) (cols x) i (rows x)))
x)
(defmethod ifft2! ((x matrix-lisp-zge))
+ (assert (and (= 1 (logcount (rows x)))
+ (= 1 (logcount (cols x)))))
+
(ifft1! x)
(dotimes (i (rows x))
(fft-radix-2-blas-complex-store! :r (matrix-store x) (cols x) i (rows x)))
Modified: src/linalg/level3-linalg-interface.lisp
==============================================================================
--- src/linalg/level3-linalg-interface.lisp (original)
+++ src/linalg/level3-linalg-interface.lisp Sun Aug 16 05:46:41 2009
@@ -35,16 +35,22 @@
(:documentation "Matrix trace (sum of diagonal elements)."))
(defgeneric mdet (matrix)
- (:documentation "Matrix determinant."))
+ (:documentation "Matrix determinant.")
+ (:method :before (m)
+ (assert (square-matrix? m))))
(defgeneric minv! (a)
(:documentation "Matrix inverse. Destructive."))
(defgeneric minv (a)
- (:documentation "Matrix inverse."))
+ (:documentation "Matrix inverse.")
+ (:method :before (m)
+ (assert (square-matrix? m))))
(defgeneric m* (a b)
- (:documentation "Matrix product."))
+ (:documentation "Matrix product.")
+ (:method :before (a b)
+ (assert (= (cols a) (rows b)))))
(defgeneric m*! (a b)
(:documentation "Matrix product. Destructive."))
Modified: src/matlisp/geev.lisp
==============================================================================
--- src/matlisp/geev.lisp (original)
+++ src/matlisp/geev.lisp Sun Aug 16 05:46:41 2009
@@ -33,14 +33,18 @@
(in-package :lisplab)
(defmethod eigenvectors ((a matrix-blas-dge))
- (destructuring-bind (evals vl-mat vr-mat)
- (dgeev (copy a) nil (mcreate a 0))
- (list evals vr-mat)))
+ (if cl-user::*lisplab-liblapack-path*
+ (destructuring-bind (evals vl-mat vr-mat)
+ (dgeev (copy a) nil (mcreate a 0))
+ (list evals vr-mat))
+ (call-next-method)))
(defmethod eigenvalues ((a matrix-blas-dge))
- (destructuring-bind (evals vl-mat vr-mat)
- (dgeev (copy a) nil nil)
- evals))
+ (if cl-user::*lisplab-liblapack-path*
+ (destructuring-bind (evals vl-mat vr-mat)
+ (dgeev (copy a) nil nil)
+ evals)
+ (call-next-method)))
(defgeneric rearrange-eigenvector-matrix (v p))
@@ -129,14 +133,18 @@
(list evec vl-mat2 vr-mat2)))))
(defmethod eigenvectors ((a matrix-zge))
- (destructuring-bind (evals vl-mat vr-mat)
- (zgeev (copy a) nil (mcreate a 0))
- (list evals vr-mat)))
+ (if cl-user::*lisplab-liblapack-path*
+ (destructuring-bind (evals vl-mat vr-mat)
+ (zgeev (copy a) nil (mcreate a 0))
+ (list evals vr-mat))
+ (call-next-method)))
(defmethod eigenvalues ((a matrix-zge))
- (destructuring-bind (evals vl-mat vr-mat)
- (zgeev (copy a) nil nil)
- evals))
+ (if cl-user::*lisplab-liblapack-path*
+ (destructuring-bind (evals vl-mat vr-mat)
+ (zgeev (copy a) nil nil)
+ evals)
+ (call-next-method)))
(defun zgeev-workspace-size (n lv? rv?)
;; Ask geev how much space it wants for the work array
Modified: src/matlisp/inv.lisp
==============================================================================
--- src/matlisp/inv.lisp (original)
+++ src/matlisp/inv.lisp Sun Aug 16 05:46:41 2009
@@ -20,41 +20,47 @@
(in-package :lisplab)
(defmethod minv! ((a matrix-blas-dge))
- (let* ((N (rows a))
- (ipiv (make-array N :element-type '(unsigned-byte 32)))
- (msg "argument A given to minv is singular to working machine precision"))
- (multiple-value-bind (_ ipiv info)
- (f77::dgetrf N N (matrix-store a) N ipiv 0)
- (declare (ignore _))
- (unless (zerop info)
- (error msg))
- (let ((work (make-array N :element-type 'double-float)))
- (multiple-value-bind (_ __ info)
- (f77::dgetri N (matrix-store a) N ipiv work N 0)
- (declare (ignore _ __))
- (unless (zerop info)
- (error msg))
- a)))))
+ (if cl-user::*lisplab-liblapack-path*
+ (let* ((N (rows a))
+ (ipiv (make-array N :element-type '(unsigned-byte 32)))
+ (msg "argument A given to minv is singular to working machine precision"))
+ (multiple-value-bind (_ ipiv info)
+ (f77::dgetrf N N (matrix-store a) N ipiv 0)
+ (declare (ignore _))
+ (unless (zerop info)
+ (error msg))
+ (let ((work (make-array N :element-type 'double-float)))
+ (multiple-value-bind (_ __ info)
+ (f77::dgetri N (matrix-store a) N ipiv work N 0)
+ (declare (ignore _ __))
+ (unless (zerop info)
+ (error msg))
+ a))))
+ ;; Othervise, call native lisp implementation
+ (call-next-method)))
(defmethod minv ((a matrix-blas-dge))
(minv! (copy a)))
(defmethod minv! ((a matrix-blas-zge))
- (let* ((N (rows a))
- (ipiv (make-array N :element-type '(unsigned-byte 32)))
- (msg "argument A given to mdiv is singular to working machine precision"))
- (multiple-value-bind (_ ipiv info)
- (f77::zgetrf N N (matrix-store a) N ipiv 0)
- (declare (ignore _))
- (unless (zerop info)
- (error msg ))
- (let ((work (make-array (* 2 N) :element-type 'double-float)))
- (multiple-value-bind (_ __ info)
- (f77::zgetri N (matrix-store a) N ipiv work N 0)
- (declare (ignore _ __))
- (unless (zerop info)
- (error msg))
- a)))))
+ (if cl-user::*lisplab-liblapack-path*
+ (let* ((N (rows a))
+ (ipiv (make-array N :element-type '(unsigned-byte 32)))
+ (msg "argument A given to mdiv is singular to working machine precision"))
+ (multiple-value-bind (_ ipiv info)
+ (f77::zgetrf N N (matrix-store a) N ipiv 0)
+ (declare (ignore _))
+ (unless (zerop info)
+ (error msg ))
+ (let ((work (make-array (* 2 N) :element-type 'double-float)))
+ (multiple-value-bind (_ __ info)
+ (f77::zgetri N (matrix-store a) N ipiv work N 0)
+ (declare (ignore _ __))
+ (unless (zerop info)
+ (error msg))
+ a))))
+ ;; Othervise, call native lisp implementation
+ (call-next-method)))
(defmethod minv ((a matrix-blas-zge))
(minv! (copy a)))
Modified: src/matlisp/mul.lisp
==============================================================================
--- src/matlisp/mul.lisp (original)
+++ src/matlisp/mul.lisp Sun Aug 16 05:46:41 2009
@@ -20,17 +20,23 @@
(in-package :lisplab)
(defmethod m* ((a matrix-blas-dge) (b matrix-blas-dge))
- (let* ((m (rows a))
- (n (cols b))
- (k (cols a))
- (c (mcreate a 0 (list m n))))
- (f77::dgemm "N" "N" m n k 1.0 (matrix-store a) m (matrix-store b) k 0.0 (matrix-store c) m)
- c))
+ (if cl-user::*lisplab-liblapack-path*
+ (let* ((m (rows a))
+ (n (cols b))
+ (k (cols a))
+ (c (mcreate a 0 (list m n))))
+ (f77::dgemm "N" "N" m n k 1.0
+ (matrix-store a) m (matrix-store b) k 0.0 (matrix-store c) m)
+ c)
+ (call-next-method)))
(defmethod m* ((a matrix-blas-zge) (b matrix-blas-zge))
- (let* ((m (rows a))
- (n (cols b))
- (k (cols a))
- (c (mcreate a 0 (list m n))))
- (f77::zgemm "N" "N" m n k #C(1.0 0.0) (matrix-store a) m (matrix-store b) k #C(0.0 0.0) (matrix-store c) m)
- c))
+ (if cl-user::*lisplab-liblapack-path*
+ (let* ((m (rows a))
+ (n (cols b))
+ (k (cols a))
+ (c (mcreate a 0 (list m n))))
+ (f77::zgemm "N" "N" m n k #C(1.0 0.0)
+ (matrix-store a) m (matrix-store b) k #C(0.0 0.0) (matrix-store c) m)
+ c)
+ (call-next-method)))
Modified: start.lisp
==============================================================================
--- start.lisp (original)
+++ start.lisp Sun Aug 16 05:46:41 2009
@@ -1,30 +1,17 @@
-;; Loads lisplab. Hack this file to fit your setup
-;; and your shared libraries.
+;; Script that loads Lisplab.
(in-package :cl-user)
-;; TODO make this part of the package
-(setf *READ-DEFAULT-FLOAT-FORMAT* 'double-float)
-
+;; Uncomment or move to .sbslrc
;; (defvar *lisplab-libblas-path* #P"/usr/lib/atlas/libblas.so.3.0")
;; (defvar *lisplab-liblapack-path* #P"/usr/lib/atlas/liblapack.so.3.0")
;; (defvar *lisplab-libfftw-path* #P"/usr/lib/libfftw3.so.3")
-(defun load-lisplab ()
- (asdf:oos 'asdf:load-op 'lisplab)
- (let ((asdf:*compile-file-failure-behaviour* :ignore))
+
+(require :lisplab)
+(let ((asdf:*compile-file-failure-behaviour* :ignore))
;; There seems to bee some compilation trouble in SBCL
;; due to type interference. Should be fixed, not just skipped.
- (asdf:oos 'asdf:load-op 'slatec))
-
- ;; The system lisplab-matlisp depends on libblas.so and liblapack.so
- ;; (asdf:oos 'asdf:load-op 'lisplab-matlisp)
-
- ;; The system lisplab-fftw depends on libfftw.so
- ;; (asdf:oos 'asdf:load-op 'lisplab-fftw)
- )
+ (require :slatec))
-;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
-(load-lisplab)
-(format t "Lisplab is loaded~%")
-;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
\ No newline at end of file
+(format t "Lisplab is loaded.~%")
More information about the lisplab-cvs
mailing list