ladata.lsp

   1 ;;; -*- mode: lisp -*-
   2 ;;; Copyright (c) 2005--2007, by A.J. Rossini <blindglobe@gmail.com>
   3 ;;; See COPYRIGHT file for any additional restrictions (BSD license).
   4 ;;; Since 1991, ANSI was finally finished.  Edited for ANSI Common Lisp.
   5
   6 ;;;; ladata -- Data handling functions for linear algebra interface
   7 ;;;;
   8 ;;;; Copyright (c) 1991, by Luke Tierney. Permission is granted for
   9 ;;;; unrestricted use.
  10
  11 ;;;
  12 ;;; Package Setup
  13 ;;;
  14
  15 ;;(in-package #:lisp-stat-basics)
  16 ;;(in-package :cl-user)
  17
  18 (defpackage :lisp-stat-linalg-data
  19   (:use :common-lisp
  20         :cffi
  21         :lisp-stat-ffi-int
  22         :lisp-stat-types
  23         :lisp-stat-sequence
  24         :lisp-stat-compound-data
  25         :lisp-stat-matrix)
  26   (:export ;; more to add
  27       +mode-in+ +mode-re+ +mode-cx+ mode-of
  28
  29       la-data-mode la-allocate la-free
  30
  31       la-get-double
  32       la-matrix la-free-matrix la-matrix-to-data la-data-to-matrix
  33       la-vector la-free-vector la-vector-to-data la-data-to-vector ))
  34
  35 (in-package :lisp-stat-linalg-data)
  36
  37
  38
  39
  40 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  41 ;;;
  42 ;;;                        Data Mode Functions
  43 ;;;
  44 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  45
  46 ;;;;
  47 ;;;; These constants need to be redefined if IN, RE or CX in linalg.h
  48 ;;;; are changed.
  49 ;;;;
  50
  51 ;;; FIXME:AJR: This is how Luke got around having appropriate
  52 ;;; approaches for Linear Algebra.  We want to cheat and instead use
  53 ;;; CLEM or MATLISP as the underlying linear algebra package.
  54
  55 (defparameter +mode-in+ 0)
  56 (defparameter +mode-re+ 1)
  57 (defparameter +mode-cx+ 2)
  58
  59 (defun mode-of (x)
  60   (etypecase x
  61              (fixnum +mode-in+)
  62              (rational +mode-re+)
  63              (float +mode-re+)
  64              (complex +mode-cx+)))
  65
  66 (defun la-data-mode (data)
  67   (let ((data (compound-data-seq data))
  68         (mode 0))
  69     (cond
  70      ((vectorp data)
  71       (let ((n (length data)))
  72         (declare (fixnum n))
  73         (dotimes (i n mode)
  74                  (declare (fixnum i))
  75                  (setf mode (max mode (mode-of (aref data i)))))))
  76      ((consp data) (dolist (x data mode) (setf mode (max mode (mode-of x)))))
  77      (t (error "bad sequence - ~s" data)))))
  78
  79
  80 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  81 ;;;;
  82 ;;;;                        Internal Allocation Funcitons
  83 ;;;;
  84 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
  85
  86 (defvar *la-allocations* nil)
  87
  88 ;;;
  89 ;;; CFFI glue for...     Storage Allocation Functions
  90 ;;;
  91
  92 (defun null-ptr-p (p) (cffi:null-pointer-p p))
  93 (defun ptr-eq (p q) (cffi:pointer-eq p q))
  94
  95 (cffi:defcfun ("la_base_allocate" ccl-la-base-allocate)
  96     :pointer (n :int) (m :int))
  97 (defun la-base-allocate (n m)
  98   (ccl-la-base-allocate n m))
  99
 100 (cffi:defcfun ("la_base_free_alloc" ccl-la-base-free-alloc)
 101     :void (p :pointer))
 102 (defun la-base-free (p)
 103   (ccl-la-base-free-alloc p))
 104
 105 (cffi:defcfun ("la_mode_size" ccl-la-mode-size)
 106     :int (x :int))
 107
 108 (defun la-mode-size (mode)
 109   (ccl-la-mode-size mode))
 110
 111 ;;;
 112 ;;;             Callbacks for Internal Storage
 113 ;;;
 114
 115 (cffi:defcallback lisp-la-allocate :void ((n :long) (m :long))
 116                   (ccl-store-ptr (la-allocate n m)))
 117 (cffi:defcfun ("register_la_allocate" register-la-allocate)
 118     :void (p :pointer))
 119 (register-la-allocate (cffi:callback lisp-la-allocate))
 120 (cffi:defcfun ("la_allocate" la)
 121     :pointer (x :int) (y :int))
 122
 123 (cffi:defcallback lisp-la-free-alloc
 124     :void ((p :pointer))
 125   (la-free p))
 126
 127 (cffi:defcfun ("register_la_free_alloc" register-la-free-alloc)
 128     :void (p :pointer))
 129 (register-la-free-alloc (cffi:callback lisp-la-free-alloc))
 130 (cffi:defcfun ("la_free_alloc" lf)
 131     :void (p :pointer))
 132
 133
 134
 135 ;;; CFFI glue for...         Storage Access Functions
 136
 137
 138 (cffi:defcfun ("la_get_integer" ccl-la-get-integer)
 139     :int (p :pointer) (i :int))
 140 (defun la-get-integer (p i)
 141   (ccl-la-get-integer p i))
 142
 143 (cffi:defcfun ("la_get_double" ccl-la-get-double)
 144     :double (p :pointer) (i :int))
 145 (defun la-get-double (p i)
 146   (ccl-la-get-double p i))
 147
 148 (cffi:defcfun ("la_get_complex_real" ccl-la-get-complex-real)
 149     :double (p :pointer) (i :int))
 150 (defun la-get-complex-real (p i)
 151   (ccl-la-get-complex-real p i))
 152
 153 (cffi:defcfun ("la_get_complex_imag" ccl-la-get-complex-imag)
 154     :double (p :pointer) (i :int))
 155 (defun la-get-complex-imag (p i)
 156   (ccl-la-get-complex-imag p i))
 157
 158 (defun la-get-complex (p i)
 159   (complex (la-get-complex-real p i) (la-get-complex-imag p i)))
 160
 161 (cffi:defcfun ("la_get_pointer" ccl-la-get-pointer)
 162     :pointer (p :pointer) (i :int))
 163 (defun la-get-pointer (p i)
 164   (ccl-la-get-pointer p i))
 165
 166
 167 ;;; CFFI glue for     Storage Mutation Functions
 168
 169 (cffi:defcfun ("la_put_integer" ccl-la-put-integer)
 170     :void (p :pointer) (i :int) (x :int))
 171 (defun la-put-integer (p i x)
 172   (ccl-la-put-integer p i x))
 173
 174 (cffi:defcfun ("la_put_double" ccl-la-put-double)
 175     :void (p :pointer) (i :int) (x :double))
 176 (defun la-put-double (p i x)
 177   (ccl-la-put-double p i (float x 1d0)))
 178
 179 (cffi:defcfun ("la_put_complex" ccl-la-put-complex)
 180     :void (p :pointer) (i :int) (x :double) (y :double))
 181 (defun la-put-complex (p i x y)
 182   (ccl-la-put-complex p i (float x 1d0) (float y 1d0)))
 183
 184 (cffi:defcfun ("la_put_pointer" ccl-la-put-pointer)
 185     :void (p :pointer) (i :int) (q :pointer))
 186 (defun la-put-pointer (p i q)
 187   (ccl-la-put-pointer p i q))
 188
 189
 190 ;; User interface (exported)
 191
 192 (defun la-allocate (n m)
 193   (let ((p (la-base-allocate n m)))
 194     (if (null-ptr-p p) (error "allocation failed"))
 195     (if (member p *la-allocations* :test #'ptr-eq)
 196       (error "pointer is already on the list"))
 197     (push p *la-allocations*)
 198     p))
 199
 200 (defun la-free (p)
 201   (when (and (not (null-ptr-p p)) (member p *la-allocations* :test #'ptr-eq))
 202         (setf *la-allocations* (delete p *la-allocations* :test #'ptr-eq))
 203         (la-base-free p)))
 204
 205 (defun la-cleanup-allocations ()
 206   (let ((allocs (copy-list *la-allocations*)))
 207     (dolist (p allocs) (la-free p))))
 208
 209 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 210 ;;;;
 211 ;;;;                       C Vector and Array Allocation
 212 ;;;;
 213 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 214
 215 (defun la-vector(n mode) (la-allocate n (la-mode-size mode)))
 216 (defun la-free-vector (v) (la-free v))
 217
 218 (defun la-matrix (n m mode)
 219   (let ((matrix (la-allocate n (la-mode-size +mode-in+))))
 220     (dotimes (i n)
 221       (la-put-pointer matrix i (la-allocate m (la-mode-size mode))))
 222     matrix))
 223
 224 (defun la-free-matrix (matrix n)
 225   (dotimes (i n)
 226     (la-free (la-get-pointer matrix i)))
 227   (la-free matrix))
 228
 229
 230 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 231 ;;;;
 232 ;;;;                   C to/from Lisp Data Conversion
 233 ;;;;
 234 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
 235
 236 (defun la-data-to-vector (data mode)
 237   (check-sequence data)
 238   (let* ((n (length data))
 239          (vec (la-vector n mode))
 240          (d (make-next-element data)))
 241     (declare (fixnum n))
 242     (cond
 243      ((= mode +mode-in+)
 244       (dotimes (i n)
 245         (declare (fixnum i))
 246         (la-put-integer vec i (get-next-element d i))))
 247      ((= mode +mode-re+)
 248       (dotimes (i n)
 249         (declare (fixnum i))
 250         (la-put-double vec i (get-next-element d i))))
 251       ((= mode +mode-cx+)
 252        (dotimes (i n)
 253          (declare (fixnum i))
 254          (let ((x (get-next-element d i)))
 255            (la-put-complex vec i (realpart x) (imagpart x))))))
 256     vec))
 257
 258 (defun la-data-to-matrix (data mode)
 259   (check-matrix data)
 260   (let* ((n (num-rows data))
 261          (m (num-cols data))
 262          (mat (la-matrix n m mode)))
 263     (declare (fixnum n m))
 264     (cond
 265      ((= mode +mode-in+)
 266       (dotimes (i n)
 267         (declare (fixnum i))
 268         (let ((vec (la-get-pointer mat i)))
 269           (dotimes (j m)
 270             (la-put-integer vec j (aref data i j))))))
 271       ((= mode +mode-re+)
 272        (dotimes (i n)
 273          (declare (fixnum i))
 274          (let ((vec (la-get-pointer mat i)))
 275            (dotimes (j m)
 276              (la-put-double vec j (aref data i j))))))
 277       ((= mode +mode-cx+)
 278        (dotimes (i n)
 279          (declare (fixnum i))
 280          (let ((vec (la-get-pointer mat i)))
 281            (dotimes (j m)
 282              (let ((x (aref data i j)))
 283                (la-put-complex vec i (realpart x) (imagpart x))))))))
 284     mat))
 285
 286 (defun la-vector-to-data (vec n mode data)
 287   (declare (fixnum n))
 288   (check-sequence data)
 289   (let ((d (make-next-element data))
 290         (gf (cond
 291               ((= mode +mode-in+) #'la-get-integer)
 292               ((= mode +mode-re+) #'la-get-double)
 293               ((= mode +mode-cx+) #'la-get-complex))))
 294     (dotimes (i n)
 295       (declare (fixnum i))
 296       (set-next-element d i (funcall gf vec i))))
 297   data)
 298
 299 (defun la-matrix-to-data (mat n m mode result)
 300   (declare (fixnum n m))
 301   (check-matrix result)
 302   (let ((gf (cond
 303               ((= mode +mode-in+) #'la-get-integer)
 304               ((= mode +mode-re+) #'la-get-double)
 305               ((= mode +mode-cx+) #'la-get-complex))))
 306     (dotimes (i n)
 307       (declare (fixnum i))
 308       (let ((vec (la-get-pointer mat i)))
 309         (dotimes (j m)
 310           (declare (fixnum j))
 311           (setf (aref result i j) (funcall gf vec j))))))
 312   result)
 313