Initial commit.

No tests, this code might not work properly.

LICENSE_1_0.txt

@@ -0,0 +1,23 @@
+Boost Software License - Version 1.0 - August 17th, 2003
+
+Permission is hereby granted, free of charge, to any person or organization
+obtaining a copy of the software and accompanying documentation covered by
+this license (the "Software") to use, reproduce, display, distribute,
+execute, and transmit the Software, and to prepare derivative works of the
+Software, and to permit third-parties to whom the Software is furnished to
+do so, all subject to the following:
+
+The copyright notices in the Software and this entire statement, including
+the above license grant, this restriction and the following disclaimer,
+must be included in all copies of the Software, in whole or in part, and
+all derivative works of the Software, unless such copies or derivative
+works are solely in the form of machine-executable object code generated by
+a source language processor.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE, TITLE AND NON-INFRINGEMENT. IN NO EVENT
+SHALL THE COPYRIGHT HOLDERS OR ANYONE DISTRIBUTING THE SOFTWARE BE LIABLE
+FOR ANY DAMAGES OR OTHER LIABILITY, WHETHER IN CONTRACT, TORT OR OTHERWISE,
+ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+DEALINGS IN THE SOFTWARE.

array-operations.asd

@@ -0,0 +1,15 @@
+;;;; array-operations.asd
+
+(asdf:defsystem #:array-operations
+  :serial t
+  :description "Simple array operations library for Common Lisp."
+  :author "Tamas K. Papp <tkpapp@gmail.com>"
+  :license "Boost Software License - Version 1.0"
+  :depends-on (#:alexandria
+               #:anaphora
+               #:let-plus)
+  :pathname #P"src/"
+  :components ((:file "package")
+               (:file "utilities")
+               (:file "displacement")
+               (:file "transformations")))

src/displacement.lisp

@@ -0,0 +1,158 @@
+;;; -*- Mode:Lisp; Syntax:ANSI-Common-Lisp; Coding:utf-8 -*-
+
+(in-package #:array-operations)
+
+;;; displacement and flattening
+
+(defun displace (array dimensions &optional (offset 0))
+  "Shorthand function for displacing an array."
+  (make-array dimensions
+              :displaced-to array
+              :displaced-index-offset offset
+              :element-type (array-element-type array)))
+
+(defun flatten (array)
+  "Return ARRAY flattened to a vector.  Will share structure."
+  (displace array (array-total-size array)))
+
+;;; subarrays
+
+(defun split (array rank)
+  "Return an array of subarrays, split off at RANK.  All subarrays are
+displaced and share structure."
+  (let ((array-rank (array-rank array)))
+    (cond
+      ((or (zerop rank) (= rank array-rank))
+       array)
+      ((< 0 rank array-rank)
+       (let* ((dimensions (array-dimensions array))
+              (result (make-array (subseq dimensions 0 rank)))
+              (sub-dimensions (subseq dimensions rank))
+              (sub-size (product sub-dimensions)))
+         (dotimes (index (array-total-size result))
+           (setf (row-major-aref result index)
+                 (displace array sub-dimensions (* index sub-size))))
+         result))
+      (t (error "Rank ~A outside [0,~A]." rank array-rank)))))
+
+(defun sub-location% (dimensions subscripts)
+  "Return (values OFFSET REMAINING-DIMENSIONS) that can be used to displace a
+row-major subarray starting at SUBSCRIPTS in an array with the given
+DIMENSIONS.  NOT EXPORTED."
+  (let+ (rev-dimensions
+         rev-subscripts
+         (tail (do ((dimensions dimensions (cdr dimensions))
+                    (subscripts subscripts (cdr subscripts)))
+                   ((not subscripts) dimensions)
+                 (assert dimensions ()
+                         "More subscripts than dimensions.")
+                 (let ((s (car subscripts))
+                       (d (car dimensions)))
+                   (declare (type fixnum d))
+                   (assert (and (integerp s) (< -1 s d)) ()
+                           "Invalid subscript.")
+                   (push s rev-subscripts)
+                   (push d rev-dimensions))))
+         (product (product tail))
+         (sum 0))
+    (declare (type fixnum product sum))
+    (mapc (lambda (d s)
+            (declare (type fixnum d s))
+            (incf sum (the fixnum (* product s)))
+            (multf product d))
+          rev-dimensions rev-subscripts)
+    (values sum tail)))
+
+(defun sub (array &rest subscripts)
+  "Given a partial list of subscripts, return the subarray that starts there,
+with all the other subscripts set to 0, dimensions inferred from the original.
+If no subscripts are given, the original array is returned.  Implemented by
+displacing, may share structure."
+  (if subscripts
+      (let+ (((&values offset dimensions)
+              (sub-location% (array-dimensions array) subscripts)))
+        (if dimensions
+            (displace array dimensions offset)
+            (apply #'aref array subscripts)))
+      array))
+
+(defun (setf sub) (value array &rest subscripts)
+  (let+ (((&values subarray atom?) (apply #'sub array subscripts)))
+    (if atom?
+        (setf (apply #'aref array subscripts) value)
+        (prog1 value
+          (assert (same-dimensions? value subarray))
+          (replace (flatten subarray) (flatten value))))))
+
+(defun partition (array start &optional (end (array-dimension array 0)))
+  "Return a subset of the array, on the first indexes between START and END."
+  (let* ((d0 (array-dimension array 0))
+         (stride (/ (array-total-size array) d0)))
+    (assert (and (<= 0 start) (< start end) (<= end d0)))
+    (displace array (cons (- end start) (cdr (array-dimensions array)))
+              (* start stride))))
+
+(defun combine (array &optional element-type)
+  "The opposite of SUBARRAYS.  If ELEMENT-TYPE is not given, it is inferred
+from the first element of array, which also determines the dimensions.  If
+that element is not an array, the original ARRAY is returned as it is."
+  (let ((first (row-major-aref array 0)))
+    (if (arrayp first)
+        (let* ((dimensions (array-dimensions array))
+               (sub-dimensions (array-dimensions first))
+               (element-type (aif element-type it (array-element-type first)))
+               (result (make-array (append dimensions sub-dimensions)
+                                   :element-type element-type))
+               (length (product dimensions))
+               (displaced (displace result (cons length sub-dimensions))))
+          (dotimes (index length)
+            (setf (sub displaced index) (row-major-aref array index)))
+          result)
+        array)))
+
+;;; subvector
+
+(defun subvec (vector start &optional (end (length vector)))
+  "Displaced vector between START and END."
+  (displace vector (- end start) start))
+
+(declaim (inline (setf subvec)))
+(defun (setf subvec) (value vector start &optional (end (length vector)))
+  ;; just a synonym for (setf subseq), defined for symmetry
+  (setf (subseq vector start end) value))
+
+;;; reshaping
+
+(defun fill-in-dimensions (dimensions size)
+  "If one of the dimensions is missing (indicated with T), replace it with a
+dimension so that the total product equals SIZE.  If that's not possible,
+signal an error.  If there are no missing dimensions, just check that the
+product equals size."
+  (let+ ((dimensions (ensure-list dimensions))
+         ((&flet missing? (dimension) (eq dimension t)))
+         missing
+         (product 1))
+    (mapc (lambda (dimension)
+            (if (missing? dimension)
+                (progn
+                  (assert (not missing) () "More than one missing dimension.")
+                  (setf missing t))
+                (progn
+                  (check-type dimension (integer 1))
+                  (multf product dimension))))
+          dimensions)
+    (if missing
+        (let+ (((&values fraction remainder) (floor size product)))
+          (assert (zerop remainder) ()
+                  "Substitution does not result in an integer.")
+          (mapcar (lambda (dimension)
+                    (if (missing? dimension) fraction dimension))
+                  dimensions))
+        dimensions)))
+
+(defun reshape (array dimensions &optional (offset 0))
+  "Reshape ARRAY using DIMENSIONS, one of which may be T which is calculated
+on demand."
+  (let* ((size (array-total-size array))
+         (dimensions (fill-in-dimensions dimensions (- size offset))))
+    (displace array dimensions offset)))

src/package.lisp

@@ -0,0 +1,29 @@
+;;;; package.lisp
+
+(defpackage #:array-operations
+  (:use #:cl #:alexandria #:anaphora #:let-plus)
+  (:nicknames #:ao)
+  (:shadow #:flatten)
+  (:export ; utilities
+   #:walk-subscripts)
+  (:export ; displacement
+   #:displace
+   #:flatten
+   #:split
+   #:sub
+   #:partition
+   #:combine
+   #:subvec
+   #:reshape)
+  (:export ; transformations
+   #:generate*
+   #:generate
+   #:permutation-repeated-index
+   #:permutation-invalid-index
+   #:valid-permutation?
+   #:complement-permutation
+   #:permute
+   #:each*
+   #:each
+   #:margin*
+   #:margin))

src/transformations.lisp

@@ -0,0 +1,105 @@
+;;; -*- Mode:Lisp; Syntax:ANSI-Common-Lisp; Coding:utf-8 -*-
+
+(in-package #:array-operations)
+
+;;; creating arrays
+
+(defun generate* (element-type function dimensions &optional arguments)
+  (aprog1 (make-array dimensions :element-type element-type)
+    (ecase arguments
+      ((nil)
+       (dotimes (position (array-total-size it))
+         (setf (row-major-aref it position)
+               (funcall function))))
+      (:position
+       (walk-subscripts (dimensions subscripts position)
+         (setf (row-major-aref it position) (funcall function position))))
+      (:subscripts
+       (walk-subscripts (dimensions subscripts position)
+         (setf (row-major-aref it position)
+               (funcall function subscripts))))
+      (:position-and-subscripts
+       (walk-subscripts (dimensions subscripts position)
+         (setf (row-major-aref it position)
+               (funcall function position subscripts)))))))
+
+(defun generate (function dimensions &optional arguments)
+  (generate* t function dimensions arguments))
+
+
+
+;;; permutations
+
+(define-condition permutation-repeated-index (error)
+  ((index :initarg :index)))
+
+(define-condition permutation-invalid-index (error)
+  ((index :initarg :index)))
+
+(defun permutation-flags% (permutation &optional (rank (length permutation)))
+  (aprog1 (make-array rank
+                      :element-type 'bit :initial-element 0)
+    (map nil (lambda (p)
+               (assert (and (integerp p) (< -1 p rank)) ()
+                       'permutation-invalid-index :index p)
+               (assert (zerop (aref it p)) ()
+                       'permutation-repeated-index :index p)
+               (setf (aref it p) 1))
+         permutation)))
+
+(defun valid-permutation? (permutation)
+  "Test if PERMUTATION is a valid permutation (of rank RANK)."
+  (every #'plusp (permutation-flags% permutation)))
+
+(defun complement-permutation (permutation rank)
+  (loop for f across (permutation-flags% permutation rank)
+        for index from 0
+        when (zerop f)
+        collect index))
+
+(defun permute (array permutation)
+  "Return an array B, where
+
+  B[b_1,...,b_n] = A[a_1,...,a_n] with a_{P[i]}=b_i
+
+A is ARRAY, and P is the axes.
+
+Permute array axes.  Elements of the sequence PERMUTATION indicate where
+that particular axis is coming from in ARRAY.  Axes in permutation can be
+repeated."
+  (assert (valid-permutation? permutation))
+  (let+ ((source-dimensions (array-dimensions array))
+         (target-dimensions (map 'list (curry #'elt source-dimensions)
+                                 permutation))
+         (target (make-array target-dimensions
+                             :element-type (array-element-type array)))
+         (buffer (make-list (array-rank array))))
+    (walk-subscripts (target-dimensions subscripts position)
+      (setf (row-major-aref target position)
+            (apply #'aref array
+                   (map-into buffer
+                             (lambda (p) (aref subscripts p))
+                             permutation))))
+    target))
+
+
+
+;;; margin
+
+(defun each* (element-type function array &rest other-arrays)
+  (aprog1 (make-array (array-dimensions array) :element-type element-type)
+    (assert (apply #'same-dimensions? array other-arrays))
+    (apply #'map-into (flatten it) function
+           (flatten array) (mapcar #'flatten other-arrays))))
+
+(defun each (function array &rest other-arrays)
+  (apply #'each* t function array other-arrays))
+
+(defun margin* (element-type function array inner
+                &optional (outer (complement-permutation inner (array-rank array))))
+  (each* element-type function
+         (split (permute array (append outer inner)) (length outer))))
+
+(defun margin (function array inner
+               &optional (outer (complement-permutation inner (array-rank array))))
+  (margin* t function array inner outer))

src/utilities.lisp

@@ -0,0 +1,47 @@
+;;; -*- Mode:Lisp; Syntax:ANSI-Common-Lisp; Coding:utf-8 -*-
+
+(in-package #:array-operations)
+
+;;; utilities used internally, not exported
+
+(defun product (dimensions)
+  "Product of elements in the argument.  NOT EXPORTED."
+  (reduce #'* dimensions))
+
+(define-modify-macro multf (&rest values) * "Multiply by the arguments")
+
+(defun same-dimensions? (array &rest arrays)
+  "Test if arguments have the same dimensions.  NOT EXPORTED."
+  (let ((dimensions (array-dimensions array)))
+    (every (lambda (array)
+             (equal dimensions (array-dimensions array)))
+           arrays)))
+
+(defmacro walk-subscripts ((dimensions subscripts
+                            &optional (position (gensym "POSITION")))
+                           &body body)
+  "Iterate over the subscripts of an array with given DIMENSIONS.  SUBSCRIPTS
+contains the current subscripts as a vector of fixnums, POSITION has the
+row-major index.  Consequences are undefined if either POSITION or SUBSCRIPTS
+is modified."
+  (check-type position symbol)
+  (check-type subscripts symbol)
+  (with-unique-names (rank last increment)
+    (once-only (dimensions)
+      `(let+ ((,rank (length ,dimensions))
+              (,dimensions (make-array ,rank
+                                       :element-type 'fixnum
+                                       :initial-contents ,dimensions))
+              (,last (1- ,rank))
+              (,subscripts (make-array ,rank
+                                       :element-type 'fixnum
+                                       :initial-element 0))
+              ((&labels ,increment (index)
+                 (unless (minusp index)
+                   (when (= (incf (aref ,subscripts index))
+                            (aref ,dimensions index))
+                     (setf (aref ,subscripts index) 0)
+                     (,increment (1- index)))))))
+         (dotimes (,position (product ,dimensions))
+           ,@body
+           (,increment ,last))))))