Skip to content
Snippets Groups Projects
Commit d0ba4920 authored by Liam Healy's avatar Liam Healy
Browse files

Default arguments for eigensystems

By default, create the appropriate output vectors and matrices, as
well as the workspace, for the eigenvalue and eigenvector functions.
parent 52686768
Branches
Tags
No related merge requests found
;; Foreign arrays (usually in C)
;; Liam Healy 2008-12-28 10:44:22EST foreign-array.lisp
;; Time-stamp: <2008-12-29 14:34:24EST foreign-array.lisp>
;; Time-stamp: <2009-01-04 21:42:09EST foreign-array.lisp>
;; $Id: $
(in-package :gsl)
......@@ -42,7 +42,8 @@
;;; Allowable keys: :dimensions, :initial-contents, :initial-element.
(defmethod initialize-instance :after
((object foreign-array) &rest initargs &key dimensions initial-contents initial-element)
((object foreign-array) &rest initargs
&key dimensions initial-contents initial-element)
(declare (ignore dimensions initial-contents initial-element))
(with-slots (cl-array dimensions original-array offset total-size) object
(let ((ffa (apply #'make-ffa (element-type object) initargs)))
......
;; Eigenvectors and eigenvalues
;; Liam Healy, Sun May 21 2006 - 19:52
;; Time-stamp: <2008-12-26 18:57:10EST eigensystems.lisp>
;; Time-stamp: <2009-01-04 21:27:10EST eigensystems.lisp>
;; $Id$
(in-package :gsl)
......@@ -43,9 +43,19 @@
;;;; Eigenvalues and eigenvectors
;;;;****************************************************************************
(defmfun eigenvalues ((A matrix) eigenvalues ws)
(defmfun eigenvalues
((A matrix)
&optional
(eigenvalues
(make-marray
(element-type A)
:dimensions (first (dimensions A))))
(ws
(if (eq (class-of A) 'complex-double-float)
(make-eigen-herm (first (dimensions A)))
(make-eigen-symm (first (dimensions A))))))
(double-float "gsl_eigen_symm"
complex-double-float "gsl_eigen_herm")
complex-double-float "gsl_eigen_herm")
(((mpointer A) :pointer)
((mpointer eigenvalues) :pointer) ((mpointer ws) :pointer))
:definition :generic
......@@ -63,9 +73,23 @@
referenced. The eigenvalues are stored in the vector eigenvalues and
are unordered.")
(defmfun eigenvalues-eigenvectors ((A matrix) eigenvalues eigenvectors ws)
(defmfun eigenvalues-eigenvectors
((A matrix)
&optional
(eigenvalues
(make-marray
(element-type A)
:dimensions (first (dimensions A))))
(eigenvectors
(make-marray
(element-type A)
:dimensions (dimensions A)))
(ws
(if (eq (class-of A) 'complex-double-float)
(make-eigen-hermv (first (dimensions A)))
(make-eigen-symmv (first (dimensions A))))))
(double-float "gsl_eigen_symmv"
complex-double-float "gsl_eigen_hermv")
complex-double-float "gsl_eigen_hermv")
(((mpointer A) :pointer) ((mpointer eigenvalues) :pointer)
((mpointer eigenvectors) :pointer) ((mpointer ws) :pointer))
:definition :generic
......@@ -118,11 +142,11 @@
(defun eigenvalue-eigenvectors-example ()
(let ((evecs (make-marray 'double-float :dimensions '(3 3)))
(evals (make-marray 'double-float :dimensions 3))
(ws (make-eigen-symmv 3))
(mat #m((20.0d0 -10.0d0 0.0d0)
(-10.0d0 30.0d0 0.0d0)
(0.0d0 0.0d0 40.0d0))))
(evals (make-marray 'double-float :dimensions 3))
(ws (make-eigen-symmv 3))
(mat #m((20.0d0 -10.0d0 0.0d0)
(-10.0d0 30.0d0 0.0d0)
(0.0d0 0.0d0 40.0d0))))
(eigenvalues-eigenvectors mat evals evecs ws)
(values (cl-array evals) (cl-array evecs))))
......
;; Definition of GSLL system
;; Liam Healy
;; Time-stamp: <2009-01-03 14:10:02EST gsll.asd>
;; Time-stamp: <2009-01-04 21:44:14EST gsll.asd>
;; $Id$
(asdf:defsystem "gsll"
......@@ -35,7 +35,7 @@
:components
((:file "foreign-friendly")
(:file "foreign-array" :depends-on (foreign-friendly))
(:file "marray" :depends-on (foreign-friendly))
(:file "marray" :depends-on (foreign-array))
(:file "vector" :depends-on (marray))
(:file "matrix" :depends-on (marray vector))
(:file "maref" :depends-on (marray vector matrix))
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment