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