Is there any good efficient vector library?

Discussion of Common Lisp

Is there any good efficient vector library?

Postby Jasper » Fri Mar 27, 2009 6:54 am

Is there any efficient vector library out there? How well is the type (vector &optional element-type size) optimized? (I don't want to sacrifice any performence.) I started making a vector package, but i don't want to duplicate effort. (Well, i had a vector library, but it doesn't satisfy points below.) I can't find any vector libraries, unless they come with other things.(So there is duplication effort in the projects there.

What i expect from a proper vector library:
  • No specifying what type of vector it is anywhere except in the type or in some cases the creation.
  • No storage of integer giving the length when lengths are small. No iterating, expand small loops.
  • v+, v-, v*, v/, inpr, crosspr, lensqr, len. angle-of-v (values if multiple.), v-of-angle (More arguments if multiple.)
  • Preferably (define-vector vector-type maker component-type initform components)-like macro that integrates with the rest. (Same rules for types as with numbers. (* float integer)->float and all. Not shy of integers.
  • Matrix stuff separate package. Might want to call matrix-multiplication the same inpr as the vectors.
Posts: 209
Joined: Fri Oct 10, 2008 8:22 am
Location: Eindhoven, The Netherlands

Re: Is there any good efficient vector library?

Postby gugamilare » Fri Mar 27, 2009 7:27 am

A quick search into gave me the following results: (for sparse vectors / matrices)

Anyway, it doesn't hurt to program your own library, just for learning and without worrying with performance, but if you want a good one, you can look into those one, and maybe a few more.

Take a look at as well.
Posts: 406
Joined: Sat Mar 07, 2009 6:17 pm
Location: Brazil

Re: Is there any good efficient vector library?

Postby Jasper » Fri Mar 27, 2009 3:28 pm

Those are a lot of libraries to check. Looking superficially none seem to have what i want.

I made a vector package. It satisfies most of my requirements. It even makes matrices as extensions of macros (inpr being matrix multiplication.) I tried to make the addition/substraction work 'in one go' rather then one by one. I mean doing (+ a b c) as such instead of (+ a (+ b c)) ) but unfortunately that made the method system go crazy, taking ages to load and such, only other way to fix it i see violates my first demand.

Code below is 318 lines total.

vect.lisp Package, basic definitions.
Code: Select all
(defpackage #:vect
  (:use #:common-lisp #:iterate)
  (:export mk-vect define-vect u v w x y z
      v+ v- v* v/ inpr crosspr lensqr len distsqr dist
      v-angle angle-v2 angle-v3
  (:documentation "Vectors/matrices with numbers in them."))

(in-package #:vect)

(defgeneric r* (a b)
  (:documentation "Multiplies two numbers. v* is the macro that users use.\
 There is no r/ it is done via r*."))

(defgeneric r+2 (a b)
  (:documentation "Adds two numbers. v+ is the macro that users use. It
leads to r+n methods."))
(defgeneric r-2 (a b)
  (:documentation "subtracts two numbers. v+ is the macro that users use.\
 It leads to r+n methods."))

(defgeneric r+3 (a b c) (:documentation "See r+2"))
(defgeneric r+4 (a b c d) (:documentation "See r+2"))
(defgeneric r+5 (a b c d e) (:documentation "See r+2"))
(defgeneric r+6 (a b c d e f) (:documentation "See r+2"))

(defgeneric r-3 (a b c) (:documentation "See r+2"))
(defgeneric r-4 (a b c d) (:documentation "See r+2"))
(defgeneric r-5 (a b c d e) (:documentation "See r+2"))
(defgeneric r-6 (a b c d e f) (:documentation "See r+2"))

(defgeneric mk-vect2 (x y)
  (:documentation "Makes two-dimensional vector, with type inferred from\

(defgeneric mk-vect3 (x y z)
  (:documentation "Makes three-dimensional vector, with type inferred from\

(defgeneric mk-vect4 (w x y z)     (:documentation "See mk-vect2"))
(defgeneric mk-vect5 (v w x y z)   (:documentation "See mk-vect2"))
(defgeneric mk-vect6 (u v w x y z) (:documentation "See mk-vect2"))

(defun vect-maker-name (n)
  (case n (2 'mk-vect2) (3 'mk-vect3) (4 'mk-vect4)
     (5 'mk-vect5) (6 'mk-vect6)))
(defmacro mk-vect (&rest args)
  `(,(vect-maker-name (length args)) ,@args))

(defmacro v* (vect &rest factors)
  "Multiplies vect with factors."
  (if (null factors) vect
      `(r* ,vect (v* ,@factors))))

(defmacro v/ (vect &rest dividers)
  "Divides vect with dividers."
  (if (null dividers) vect
      `(r* ,vect (/ (v* ,@dividers)))))

(defmacro v+ (&rest add)
  "Addition macro. Chooses from the r+n functions, chains them if more then\
  (case (length add)
    (1 (car add))    (2 `(r+2 ,@add)); (3 `(r+3 ,@add))
;    (4 `(r+4 ,@add)) (5 `(r+5 ,@add)) (6 `(r+6 ,@add))
    (t `(v+ (v+2 ,@(subseq add 0 2)) ,@(subseq add 2)))))

(defmacro v- (&rest add)
  "Substraction macro. Chooses from the r-n functions, chains them if\
 more then 6."
  (case (length add)
    (1 (car add))    (2 `(r-2 ,@add)); (3 `(r-3 ,@add))
    (t `(v- (v-2 ,@(subseq add 0 2)) ,@(subseq add 2)))))

;    (4 `(r-4 ,@add)) (5 `(r-5 ,@add)) (6 `(r-6 ,@add))
;    (t `(v- (v-6 ,@(subseq add 0 6)) ,@(subseq add 6)))))

(defgeneric r-inpr (a b) (:documentation "Inproduct between two things.
Used through inpr macro."))

(defgeneric crosspr (a b) (:documentation "Cross product between two\

(defmacro inpr (&rest args)
  "Inproduct between things. Usually two."
  `(r-inpr ,(car args) (inpr ,@(cdr args))))

(declaim (inline len lensqr distsqr dist))
(defun lensqr (vect)
  "Length squared of a vector."
  (r-inpr vect vect))

(defun len (vect)
  "Length of a vector."
  (sqrt (lensqr vect)))

(defun distsqr (a b)
  "Distance squared between two vectors."
  (lensqr (v- a b)))
(defun dist (a b)
  "Distance between two vectors."
  (sqrt (distsqr a b)))

(defgeneric v-angle (vect) (:documentation "Angles from a vector."))

(defun angle-v2 (angle)
  "Vector from a angle."
  (mk-vect (cos angle) (sin angle)))

(defun angle-v3 (angle up-angle)
  "Vector from two angles, so three-dimensional.
 (Not the usual illogical polar coordinates, height is sin up-angle.)"
  (mk-vect (* (cos angle) (cos up-angle)) (* (sin angle) (cos up-angle))
      (sin up-angle)))

(defgeneric vect-as-list (vec) (:documentation "Converts a vector to a list\
 for printing and such."))
vect-number.lisp The vector code treats things internally as 'vectors' so that it can do matrices as well, hence we need to make a few methods to make numbers work.
Code: Select all
(in-package #:vect)

;;Regular number version.
(defmethod r* ((x number) (y number))
  (* x y))

;Addition and substraction.
(eval (let*((vars '(v1 v2 v3 v4 v5 v6))
        (args (iter (for v in vars)
          (collect `(,v number)))))
    `(progn ,@(iter (for n from 2)
          (for r+ in '(r+2 r+3 r+4 r+5 r+6))
          (for r- in '(r-2 r-3 r-4 r-5 r-6))
          (appending `((defmethod ,r+ (,@(subseq args 0 n))
               (+ ,@(subseq vars 0 n)))
                  (defmethod ,r- (,@(subseq args 0 n))
               (+ ,@(subseq vars 0 n)))))))))
(defmethod r-inpr ((a number) (b number))
  (* a b))
vect-define.lisp Allows one to make vectors, it connects it to other vectors via mk-vect (which can be overridden)
Code: Select all
(in-package #:vect)

(defun components-from-integer (dimension)
  "Uses standard component names."
  (case dimension (1 '(x)) (2 '(x y)) (3 '(x y z))))

(defmacro do-components (fun maker components &rest args)
  "Macro, does function with the arguments from all the arguments.
Components is the number of elements, or their names."
  (when (integerp components)
    (setf components (components-from-integer)))
    ,@(iter (for c in components)
        ,@(iter (for a in args)
          (collect `(,c ,a))))))))

(defmacro make-do-components-method (method fun maker components type
  "Makes a method from the components."
  (let ((gs (iter (repeat arg-cnt)
        (collect (gensym)))))
    `(defmethod ,method (,@(iter (for g in gs)
             (collect `(,g ,type))))
       (do-components ,fun ,maker ,components ,@gs))))

(defun any-eql (eql list)
  "Any of list eql to eql."
  (dolist (el list) (when (eql eql el) (return el))))

(defun all-combinations (list n fun &optional args)
  (dolist (el list)
    (if (= n 1)
      (funcall fun (cons el args))
      (all-combinations list (- n 1) fun (cons el args)))))

(defun r+n-name (n) (case n (2 'r+2) (3 'r+3) (4 'r+4) (5 'r+5) (6 'r+6)))
(defun r-n-name (n) (case n (2 'r-2) (3 'r-3) (4 'r-4) (5 'r-5) (6 'r-6)))

(defvar *vector-types* nil "Types that already exist.")

(defmacro define-vect (name of-type initform components &key maker
  "Defines vector type with components of given type.
Also makes all methods."
  (when (integerp components)
    (setf components (components-from-integer components)))
  (unless maker
    (setf maker (vect-maker-name (length components))))
  (when (any-eql name *vector-types*)
    (error "Vector type with this name already exists."))

  (push name *vector-types*)
   ;The class.
     (defclass ,name ()
       (,@(iter (for c in components)
      (collect `(,c :initarg ,c :accessor ,c
               :type ,of-type :initform ,initform))))
       (:documentation ,(format nil "A ~D-dimensional vector of type ~D."
            (length components) of-type)))
     (defmethod ,maker (,@(iter (for c in components)
            (collect `(,c ,of-type))))
       (make-instance ',name ,@(iter (for c in components)
                 (appending `(',c ,c)))))
      ,@(when fun-maker
    `((defun ,fun-maker (,@components)
        (make-instance ',name ,@(iter (for c in components)
                  (appending `(',c ,c)))))))
   ;Multiplication. (Dividing based on it.)
     (defmethod r* ((vect ,name) (scalar number))
       (mk-vect ,@(iter (for c in components)
         (collect `(* (,c vect) scalar)))))
   ;Adding/substracting, different counts of arguments.
     ,@(iter (for n from 2) (while (<= n 2))
        (all-combinations *vector-types* n
       (lambda (tps)
         (when (any-eql name tps)
           (let*((vars (iter (repeat n) (collect (gensym))))
            (args (iter (for v in vars)
                   (for tp in tps)
                   (collect `(,v ,tp)))))
         `((defmethod ,(r+n-name n) (,@args)
             (do-components v+ mk-vect ,components
           (defmethod ,(r-n-name n) (,@args)
             (do-components v- mk-vect ,components
   ;Other stuff.
   (all-combinations *vector-types* 2
     (lambda (tps)
       (when (any-eql name tps)
         (destructuring-bind (tp-a tp-b) tps
       `((defmethod r-inpr ((a ,tp-a) (b ,tp-b)) ;Inproduct
           (v+ ,@(iter (for c in components)
             (collect `(r-inpr (,c a) (,c b))))))
         ,(case (length components)
            (2 ;Cross product of 2 dimensional is 1 dimensional.
             (destructuring-bind (x y) components
          `(defmethod crosspr  ((a ,tp-a) (b ,tp-b))
             (- (* (,y a) (,x b)) (* (,x a) (,y b))))))
             (destructuring-bind (x y z) components
          `(defmethod crosspr  ((a ,tp-a) (b ,tp-b))
              (- (* (,z a) (,y b)) (* (,y a) (,z b)))
              (- (* (,x a) (,z b)) (* (,z a) (,x b)))
              (- (* (,y a) (,x b)) (* (,x a) (,y b))))))))))))))
     (defmethod vect-as-list ((vec ,name))
       ,(if (any-eql of-type (cdr *vector-types*))
       `(iter (for el in (do-components identity list ,components vec))
         (collect (vect-as-list el)))
       `(do-components identity list ,components vec)))
     ;TODO v-angle
vect-specific.lisp Some vectors themselves.
Code: Select all
(in-package #:vect)

;;Matrices from vectors.
(define-vect vect2d double-float 0d0 2 :fun-maker mk-vect2d)
(define-vect vect2f single-float 0.0 2 :fun-maker mk-vect2f)
(define-vect vect2i fixnum       0   2 :fun-maker mk-vect2i)

(define-vect vect3d double-float 0d0 3 :fun-maker mk-vect3d)
(define-vect vect3f single-float 0.0 3 :fun-maker mk-vect3f)
(define-vect vect3i fixnum       0   3 :fun-maker mk-vect3i)
Some matrices made from those
Code: Select all
(in-package #:vect)

;;Matrices from vectors.
(define-vect matrix22d vect2d (mk-vect2d 0d0 0d0) 2
        :fun-maker mk-matrix22d)
(define-vect matrix23d vect3d (mk-vect3d 0d0 0d0 0d0) 2
          :fun-maker mk-matrix23d)
(define-vect matrix32d vect2d (mk-vect2d 0d0 0d0) 3
          :fun-maker mk-matrix32d)
(define-vect matrix33d vect3d (mk-vect3d 0d0 0d0 0d0) 3
          :fun-maker mk-matrix33d)
Posts: 209
Joined: Fri Oct 10, 2008 8:22 am
Location: Eindhoven, The Netherlands

Re: Is there any good efficient vector library?

Postby Kohath » Sat Mar 27, 2010 4:14 pm

Hey, as a note (and I know it's exactly one year from the last post :mrgreen:, so things may have changed since then), I would suggest implementing the r+ and r- functionality with something like this:
Code: Select all
(defun v+ (&rest add)
  (reduce #'r+2 (cdr add) :initial-value (car add)))

(defun v- (&rest add)
  (cond ((cdr add) (reduce #'r-2 (cdr add) :initial-value (car add)))
        (add (negative-of (car add)))
        (t default-value-or-error)))

The function v+ handles the no arguments case by returning nil, but an (or (car add) default-value-or-error) would be one way to fix that if it's a problem.

But if you wanted to do it in one swoop, without calling the 2 argument version, I'd recommend an iterative version, checking the types of each argument and then only if all types are the same (or compatible) would I dispatch to a generic function (which would use the type of one arg to specialize), like this:
Code: Select all
(defmethod same-type-vmult ((prototype (vector integer)) vector-list)
  (loop for vect in vector-list
    do ....))
User avatar
Posts: 61
Joined: Mon Jul 07, 2008 8:06 pm
Location: Toowoomba, Queensland, Australia

Re: Is there any good efficient vector library?

Postby Jasper » Mon Mar 29, 2010 6:53 am

The newer code can be found in the vect/ directory of the git of this project. Note that i haven't done much at all to that project lately! Infact, i moved some of it to the Gil project so people wouldn't have to download both, and because i can really consider :generic so a lib for other people to use. Some of it is more stuff i like/use and will hopefully come out of a nicely distributed package somewhere later. Consider that lisp-umac thing dead. (Perhaps i should put that on the website)

I would write the code in vect/ rather differently now, although i don't see much way around it. Imo whether it is better to use classes or cl:vectors is pretty much up to the question 'which is optimized better', these days i use vector-vect and hope it will be optimized.

I think having a library(/asd package) focussing solely vanilla vector/matrix operations would be a good idea. vector-vect.lisp in that git repo doesn't do matrices and vectors in one sweep; (for convenience)i'll just copy code that does that below.
Code: Select all
;;Placed in public domain.

(cl:in-package :cl-user)

(defpackage :vector-vect
  (:use :common-lisp)
  (:export vect2d vect3d x y z xyz xyzs with-xyz with-xyzs
      v+ v- v* v/
      inpr crosspr lensqr len normalize  distsqr dist
      angle-v v2-angle rotate-v2 rotate-v3)
  (:documentation "Vectors from sequences.
Hope optimization will do its job well!"))

(in-package :vector-vect)

(deftype vect2d () '(vector number 2))
(deftype vect3d () '(vector number 3))

(defmacro mk-component (n name)
  "Simple utility to make the component functions."
  `(defmacro ,name (vect)
     ,(format nil "~D component of vector. (index ~D)" name n)
     `(aref ,vect ,',n)))

(mk-component 0 x)
(mk-component 1 y)
(mk-component 2 z)

(defmacro with-gensyms ((&rest vars)&body body)
  "Makes you some variables with gensyms output in them."
  `(let (,@(loop for el in vars
         collect `(,el (gensym))))

(defmacro xyz (input &body body)
  "Makes settable x,y,z for you, or name-x,y,z"
  (with-gensyms (v)
    (multiple-value-bind (vect name)
   (typecase input
     (list (values (car input) (cadr input)))
     (t    (values input nil)))
      (typecase name
   (null ;No name
    `(let ((,v ,vect))
      ((x (aref ,v 0)) (y (aref ,v 1)) (z (aref ,v 2)))
 ;Automatic name. (Note: not sure if providing it is good idea.)
   ((eql t)
    (typecase vect
      ((symbolp vect) ;Redo it with variable name as symbol.
       `(xyz (,vect ,vect) ,@body))
      ((listp vect) ;Use function/macro name to identify.
       `(xyx (,vect ,(car vect)) ,@body))
       (error "Don't know how to automatically give a name"))))
   ((or symbol string) ;Given name.
         (,@(loop for c in '(x y z)
          for k from 0
          collect `(,(intern(format nil "~D-~D" name c))
                (aref ,vect ,k))))
    (error "Didn't know what to do with input."))))))

(defmacro xyzs (inputs &body body)
  "Multiple xyz macros in sequence."
  (if (null (cdr inputs))
    `(xyz ,(car inputs) ,@body)
    `(xyz ,(car inputs) (multi-xyz ,(cdr inputs) ,@body))))

(defmacro with-xyz (input &body body)
  "Synonym for xyz"
  `(xyz ,input ,@body))
(defmacro with-xyzs (inputs &body body)
  "Synonym for xyzs"
  `(xyzs ,inputs ,@body))

(defmacro with-make-array (var dim k &body body)
  `(let ((,var (make-array ,dim :initial-element 0)))
     (dotimes (,k ,dim) ,@body) ;Note abstraction leak on dim.

;;Adding, substracting.
(defun v+ (&rest vects)
  (reduce (lambda (a b)
       (typecase a
          (with-make-array out (length a) i
       (setf (aref out i) (v+ (aref a i) (aref b i)))))
          (+ a b))))

(defun v* (vect &rest factors)
  "Multiplication by _scalars_"
  (let ((factor (reduce #'* factors)))
    (typecase vect
       (with-make-array out (length vect) i
    (setf (aref out i) (* (aref vect i) factor))))
       (* vect factor)))))

(defun v- (&rest vects)
  (v+ (car vects)
      (v* (apply #'v+ (cdr vects)) -1)))

(defun v/ (vect &rest divide-by)
  "Division by scalars."
  (v* vect (/ (reduce divide-by))))

(declaim (inline v+ v- v* v/))

(defun inpr (&rest vals)
  "Inproducts and matrix multiplication."
  (reduce (lambda (a b)
       (typecase a
          (print (list a b))
          (typecase b
        (do*((i 0 (+ i 1))
             (out (inpr (aref a 0) (aref b 0))
             (v+ out (inpr (aref a i) (aref b i)))))
            ((>= (+ i 1) (length a)) out)))
        (v* a b))))
          (* a b))))

(defun crosspr (a b)
  "Cross product. Only 2 or 3 dimensional."
  (xyz (a)
    (case (length a)
      (2 (- (* x (y b)) (* y (x b))))
      (3 (vector
     (- (* y (z b)) (* z (y b)))
     (- (* z (x b)) (* x (z b)))
     (- (* x (y b)) (* y (x b))))))))

(defun lensqr (&rest vects)
  "Length squared of vector."
  (loop for v in vects sum (inpr v v)))
(defun len (&rest vects)
  "Length of vector."
  (sqrt (apply #'lensqr v)))

(defun normalize (v)
  (v/ v (len v)))

(declaim (inline inpr crosspr lensqr len normalize))

;;A little bit arbitrary to have this one generic.. the others can be seen abstractly too..
(defgeneric distsqr (a b)
  (:documentation "Distance squared between two things."))
(defmethod distsqr ((a vector) (b vector))
  (lensqr (v- a b)))
(defgeneric dist (a b)
  (:documentation "Distance between two things."))
(defmethod dist (a b)
  (sqrt (distsqr a b)))

(defun angle-v-raw (angles &key pre)
  (if (null angles)
    `((* ,@pre))
    `((* ,@pre (cos ,(car angles)))
      ,@(angle-v-raw (cdr angles) :pre `(,@pre (sin ,(car angles)))))))

(defmacro angle-v (&rest angles)
  "Produces an n+1 dimensional vector from n angles."
  (let ((angle-vars (loop for el in angles
          collect (gensym))))
    `(let (,@(loop for el in angles
         for gs in angle-vars
      collect `(,gs ,el)))
       `(vector ,@(angle-v-raw angle-vars)))))

(defun v2-angle (v)
  "Gets angle from 2d vector."
  (let ((len (len v)))
    (xyz v
      (if (= len 0) 0
   (* (acos (/ x len)) (if (< y 0) -1 1))))))

(defun rotate-v2 (v a &key (cos-a (cos a)) (sin-a (sin a)))
  "Rotates a 2d vector."
  (xyz v
    (vector (+ (* y sin-a) (* x cos-a))
       (- (* y cos-a) (* x sin-a)))))

(defun rotate-v3 (v n a &key (cos-a (cos a)) (sin-a (sin a)))
  "Rotate 3d vector round something"
  (declare (type vect3d v n)
      (type number a cos-a sin-a))
  (let*((e1 (xyz (n)
       ((and (= x 0) (= y 0))
        (return-from rotate-v3
          (xyz (v)
            (vector (+ (* y sin-a) (* x cos-a))
               (- (* y cos-a) (* x sin-a))
       ((= x 0)
        (vector 0d0 (/ -1 y) (/ z)))
       ((= y 0)
        (vector (/ -1 x) 0d0 (/ z)))
        (vector (/ x) (/ -1 y) 0d0))))))
   (e2 (normalize (crosspr e1 n)))
   (x (inpr e1 v))
   (y (inpr e2 v)))
    (declare (type vect3d e1 e2)
        (type number x y))
    (v+ (v* e1 (+ (* y sin-a) (* x cos-a)))
   (v* e2 (- (* y cos-a) (* x sin-a)))
   (v* n  (inpr (normalize n) v)))))

(defun rotate-matrix-2 (a &key (cos-a (cos a)) (sin-a (sin a)))
  "Make the matrix to rotate"
  (declare (type number a cos-a cos-b))
  (vector (vector cos-a     sin-a)
     (vector (- sin-a) cos-a)))

;TODO defun rotate-matrix-3, matrix stuff etcetera.
Posts: 209
Joined: Fri Oct 10, 2008 8:22 am
Location: Eindhoven, The Netherlands

Return to Common Lisp

Who is online

Users browsing this forum: No registered users and 12 guests