Fcalva/chibi-scheme
1
0
Fork 0
mirror of https://github.com/ashinn/chibi-scheme.git synced 2025-07-04 11:46:37 +02:00
Code Issues Projects Releases Packages Wiki Activity Actions

Initial SRFI 231 implementation.

Browse source
This commit is contained in:
Alex Shinn 2022-06-12 16:04:38 +09:00
parent 4d0ae090b7
commit 38fc7e0932
5 changed files with 5088 additions and 0 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

53
lib/srfi/231.sld Normal file
View file

@ -0,0 +1,53 @@
(define-library (srfi 231)
(import (scheme base)
(scheme list)
(scheme vector)
(scheme sort)
(srfi 160 base)
(srfi 231 base)
(chibi assert))
(export
;; Miscellaneous Functions
translation? permutation?
;; Indexes
index-rotate index-first index-last
;; Intervals
make-interval interval? interval-dimension interval-lower-bound
interval-upper-bound interval-lower-bounds->list
interval-upper-bounds->list interval-lower-bounds->vector
interval-upper-bounds->vector interval= interval-volume
interval-subset? interval-contains-multi-index? interval-projections
interval-for-each interval-dilate interval-intersect
interval-translate interval-permute
interval-scale interval-cartesian-product
interval-width interval-widths
interval-empty?
;; Storage Classes
make-storage-class storage-class? storage-class-getter
storage-class-setter storage-class-checker storage-class-maker
storage-class-copier storage-class-length storage-class-default
generic-storage-class s8-storage-class s16-storage-class
s32-storage-class s64-storage-class u1-storage-class
u8-storage-class u16-storage-class u32-storage-class
u64-storage-class f8-storage-class f16-storage-class
f32-storage-class f64-storage-class
c64-storage-class c128-storage-class
storage-class-data? storage-class-data->body
;; Arrays
make-array array? array-domain array-getter array-dimension
mutable-array? array-setter specialized-array-default-safe?
specialized-array-default-mutable? make-specialized-array
make-specialized-array-from-data
specialized-array? array-storage-class array-indexer array-body
array-safe? array-packed? specialized-array-share
array-copy array-curry array-extract array-tile array-translate
array-permute array-reverse array-sample
array-outer-product array-map array-for-each array-foldl
array-foldr array-reduce array-any array-every
array-inner-product array-stack array-append array-block
array->list list->array array->vector vector->array
array->list* list*->array array->vector* vector*->array
array-assign! array-ref array-set! array-decurry
specialized-array-reshape
)
(include "231/transforms.scm"))

496
lib/srfi/231/base.scm Normal file
View file

@ -0,0 +1,496 @@
;; Miscellaneous Functions
(define (translation? x)
(and (vector? x) (vector-every exact-integer? x)))
(define (permutation? x)
(and (translation? x)
(let* ((len (vector-length x))
(seen (make-vector len 0)))
(let lp ((i 0))
(or (>= i len)
(and (< -1 (vector-ref x i) len)
(zero? (vector-ref seen (vector-ref x i)))
(begin
(vector-set! seen (vector-ref x i) 1)
(lp (+ i 1)))))))))
(define (index-rotate n k)
(list->vector (append (iota (- n k) k) (iota k))))
(define (index-first n k)
(list->vector (cons k
(append (iota k)
(iota (- n (+ k 1)) (+ k 1))))))
(define (index-last n k)
(list->vector (append (iota k)
(iota (- n (+ k 1)) (+ k 1))
(list k))))
;; Intervals
(define-record-type Interval
(%%make-interval lb ub)
interval?
(lb interval-lb)
(ub interval-ub))
(define (%make-interval lo hi)
(assert (and (translation? lo)
(translation? hi)
(= (vector-length lo) (vector-length hi))
(vector-every <= lo hi)))
(%%make-interval lo hi))
(define (make-interval x . o)
(if (pair? o)
(%make-interval x (car o))
(%make-interval (make-vector (vector-length x) 0) x)))
(define (interval-dimension iv)
(vector-length (interval-lb iv)))
(define (interval-lower-bound iv i) (vector-ref (interval-lb iv) i))
(define (interval-upper-bound iv i) (vector-ref (interval-ub iv) i))
(define (interval-lower-bounds->list iv) (vector->list (interval-lb iv)))
(define (interval-upper-bounds->list iv) (vector->list (interval-ub iv)))
(define (interval-lower-bounds->vector iv) (vector-copy (interval-lb iv)))
(define (interval-upper-bounds->vector iv) (vector-copy (interval-ub iv)))
(define (interval-width iv i)
(- (interval-upper-bound iv i) (interval-lower-bound iv i)))
(define (interval-widths iv)
(vector-map - (interval-lb iv) (interval-ub iv)))
(define (interval= iv1 iv2)
(assert (and (interval? iv1) (interval? iv2)))
(and (equal? (interval-lb iv1) (interval-lb iv2))
(equal? (interval-ub iv1) (interval-ub iv2))))
(define (interval-volume iv)
(vector-fold (lambda (acc lower upper) (* acc (- upper lower)))
1
(interval-lb iv) (interval-ub iv)))
(define (interval-empty? iv)
(zero? (interval-volume iv)))
(define (interval-subset? iv1 iv2)
(assert (and (interval? iv1) (interval? iv2)
(= (interval-dimension iv1) (interval-dimension iv2))))
(and (vector-every >= (interval-lb iv1) (interval-lb iv2))
(vector-every <= (interval-ub iv1) (interval-ub iv2))))
(define (interval-contains-multi-index? iv i0 . o)
(assert (interval? iv))
(let ((i (list->vector (cons i0 o))))
(assert (and (= (interval-dimension iv) (vector-length i))
(vector-every integer? i)))
(and (vector-every >= i (interval-lb iv))
(vector-every < i (interval-ub iv)))))
(define (interval-projections iv rd)
(values (make-interval (vector-copy (interval-lb iv) 0 rd)
(vector-copy (interval-ub iv) 0 rd))
(make-interval (vector-copy (interval-lb iv) rd)
(vector-copy (interval-ub iv) rd))))
(define (rev-index-next! rev-index rev-lowers rev-uppers)
(cond
((null? rev-index) #f)
((< (caar rev-index) (- (car rev-uppers) 1))
(set-car! (car rev-index) (+ 1 (caar rev-index)))
#t)
(else
(set-car! (car rev-index) (car rev-lowers))
(rev-index-next! (cdr rev-index) (cdr rev-lowers) (cdr rev-uppers)))))
(define (interval-cursor iv)
(let* ((rev-lowers (reverse (interval-lower-bounds->list iv)))
(rev-uppers (reverse (interval-upper-bounds->list iv)))
(multi-index (interval-lower-bounds->list iv))
(rev-index (pair-fold cons '() multi-index)))
(vector multi-index rev-index rev-lowers rev-uppers)))
(define (interval-cursor-get ivc)
(vector-ref ivc 0))
(define (interval-cursor-next! ivc)
(and (rev-index-next! (vector-ref ivc 1)
(vector-ref ivc 2)
(vector-ref ivc 3))
(vector-ref ivc 0)))
(define (interval-cursor-next ivc)
(let* ((multi-index (list-copy (vector-ref ivc 0)))
(ivc (vector multi-index
(pair-fold cons '() multi-index)
(vector-ref ivc 2)
(vector-ref ivc 3))))
(and (rev-index-next! (vector-ref ivc 1)
(vector-ref ivc 2)
(vector-ref ivc 3))
(values ivc (vector-ref ivc 0)))))
(define (interval-fold kons knil iv)
(case (interval-dimension iv)
((1)
(let ((end (interval-upper-bound iv 0)))
(do ((i (interval-lower-bound iv 0) (+ i 1))
(acc knil (kons acc i)))
((>= i end) acc))))
((2)
(let ((end0 (interval-upper-bound iv 0))
(start1 (interval-lower-bound iv 1))
(end1 (interval-upper-bound iv 1)))
(do ((i (interval-lower-bound iv 0) (+ i 1))
(acc knil
(do ((j start1 (+ j 1))
(acc acc (kons acc i j)))
((>= j end1) acc))))
((>= i end0) acc))))
(else
(let ((ivc (interval-cursor iv)))
(let lp ((acc knil))
(let ((acc (apply kons acc (interval-cursor-get ivc))))
(if (interval-cursor-next! ivc)
(lp acc)
acc)))))))
(define (interval-for-each f iv)
(interval-fold (lambda (acc . multi-index) (apply f multi-index)) #f iv)
(if #f #f))
(define (interval-dilate iv lower-diffs upper-diffs)
(assert (= (interval-dimension iv)
(vector-length lower-diffs)
(vector-length upper-diffs)))
(make-interval (vector-map + (interval-lb iv) lower-diffs)
(vector-map + (interval-ub iv) upper-diffs)))
(define (interval-intersect iv0 . o)
(let ((ls (cons iv0 o)))
(assert (and (every interval? ls)
(or (null? o) (apply = (map interval-dimension ls)))))
(let ((lower (apply vector-map max (map interval-lb ls)))
(upper (apply vector-map min (map interval-ub ls))))
(and (vector-every < lower upper)
(make-interval lower upper)))))
(define (interval-translate iv translation)
(assert (translation? translation))
(interval-dilate iv translation translation))
(define (interval-permute iv perm)
(assert (and (interval? iv) (permutation? perm)))
(let* ((len (interval-dimension iv))
(lower (make-vector len))
(upper (make-vector len)))
(assert (= len (vector-length perm)))
(do ((i 0 (+ i 1)))
((>= i len) (make-interval lower upper))
(vector-set! lower i (interval-lower-bound iv (vector-ref perm i)))
(vector-set! upper i (interval-upper-bound iv (vector-ref perm i))))))
(define (interval-scale iv scales)
(assert (and (interval? iv)
(vector? scales)
(= (interval-dimension iv) (vector-length scales))
(vector-every exact-integer? scales)
(vector-every positive? scales)))
(make-interval
(vector-map (lambda (u s) (exact (ceiling (/ u s))))
(interval-ub iv)
scales)))
(define (interval-cartesian-product iv0 . o)
(make-interval (apply vector-append (map interval-lb (cons iv0 o)))
(apply vector-append (map interval-ub (cons iv0 o)))))
;; Storage Classes
(define-record-type Storage-Class
(make-storage-class getter setter checker maker copier length default data? data->body)
storage-class?
(getter storage-class-getter)
(setter storage-class-setter)
(checker storage-class-checker)
(maker storage-class-maker)
(copier storage-class-copier)
(length storage-class-length)
(default storage-class-default)
(data? storage-class-data?)
(data->body storage-class-data->body))
(define generic-storage-class
(make-storage-class
vector-ref vector-set! (lambda (x) #t) make-vector
vector-copy! vector-length #f (lambda (data) #t) (lambda (data) data)))
;; Parameters
;; Note safety is ignored in this implementation.
(define specialized-array-default-safe?
(make-parameter #f (lambda (x) (assert (boolean? x)) x)))
(define specialized-array-default-mutable?
(make-parameter #t (lambda (x) (assert (boolean? x)) x)))
;; Arrays
(define-record-type Array
(%%make-array domain getter setter storage body coeffs indexer safe? adjacent?)
array?
(domain array-domain)
(getter array-getter)
(setter array-setter %array-setter-set!)
(storage array-storage-class)
(body array-body)
(coeffs array-coeffs)
(indexer array-indexer)
(safe? array-safe?)
(adjacent? array-adjacent? array-adjacent?-set!))
(define (%make-array domain getter setter storage body coeffs
indexer safe? adjacent?)
(assert (and (interval? domain)
(procedure? getter)
(or (not setter) (procedure? setter))
(or (not storage) (storage-class? storage))))
(%%make-array
domain getter setter storage body coeffs indexer safe? adjacent?))
(define (make-array domain getter . o)
(assert (and (interval? domain) (procedure? getter)))
(%make-array domain getter (and (pair? o) (car o)) #f #f #f #f #f #f))
(define (array-dimension a)
(interval-dimension (array-domain a)))
(define (mutable-array? x)
(and (array? x) (array-setter x) #t))
(define (array-ref array . multi-index)
(apply (array-getter array) multi-index))
(define (array-set! array val . multi-index)
(apply (array-setter array) val multi-index))
(define (specialized-getter body indexer getter)
(lambda multi-index
(getter body (apply indexer multi-index))))
(define (specialized-setter body indexer setter)
(lambda (val . multi-index)
(setter body (apply indexer multi-index) val)))
;; Indexing
(define (indexer->coeffs indexer domain . o)
(let* ((verify? (and (pair? o) (car o)))
(res (make-vector (+ 1 (interval-dimension domain)) 0))
(multi-index (interval-lower-bounds->list domain))
(base (apply indexer multi-index)))
(vector-set! res 0 base)
(let lp ((i 1)
(ls multi-index)
(offset base)
(count 0))
(cond
((null? ls)
(if (and verify? (zero? count))
(lp 1 multi-index offset (+ count 1))
res))
((= (+ 1 (interval-lower-bound domain (- i 1)))
(interval-upper-bound domain (- i 1)))
(lp (+ i 1) (cdr ls) offset count))
(else
(let ((dir (if (and (> count 0)
(= (+ (car ls) 1)
(interval-upper-bound domain (- i 1))))
-1
1)))
(set-car! ls (+ (car ls) dir))
(let* ((offset2 (apply indexer multi-index))
(coeff (* dir (- offset2 offset))))
(cond
((> count 0)
(and (= coeff (vector-ref res i))
(lp (+ i 1) (cdr ls) offset2 count)))
(else
(vector-set! res i coeff)
(vector-set! res 0 (- (vector-ref res 0)
(* coeff
(interval-lower-bound domain (- i 1)))))
(lp (+ i 1) (cdr ls) offset2 count))))))))))
(define (coeffs->indexer coeffs domain)
(case (vector-length coeffs)
((2)
(let ((a (vector-ref coeffs 0))
(b (vector-ref coeffs 1)))
(lambda (x) (+ a (* b x)))))
((3)
(let ((a (vector-ref coeffs 0))
(b (vector-ref coeffs 1))
(c (vector-ref coeffs 2)))
(lambda (x y) (+ a (* b x) (* c y)))))
((4)
(let ((a (vector-ref coeffs 0))
(b (vector-ref coeffs 1))
(c (vector-ref coeffs 2))
(d (vector-ref coeffs 3)))
(lambda (x y z) (+ a (* b x) (* c y) (* d z)))))
(else
(lambda multi-index
(let ((lim (vector-length coeffs)))
(let lp ((ls multi-index)
(i 1)
(res (vector-ref coeffs 0)))
(cond
((null? ls)
(if (< i lim)
(error "multi-index too short for domain" multi-index domain)
res))
((>= i lim)
(error "multi-index too long for domain" multi-index domain))
(else
(lp (cdr ls)
(+ i 1)
(+ res (* (car ls) (vector-ref coeffs i))))))))))))
(define (default-coeffs domain)
(let* ((dim (interval-dimension domain))
(res (make-vector (+ 1 dim))))
(vector-set! res 0 0)
(vector-set! res dim 1)
(let lp ((i (- dim 1))
(scale 1))
(cond
((< i 0)
res)
((= (+ 1 (interval-lower-bound domain i))
(interval-upper-bound domain i))
(vector-set! res (+ i 1) 0)
(lp (- i 1) scale))
(else
(let ((coeff (* scale (- (interval-upper-bound domain i)
(interval-lower-bound domain i)))))
(vector-set! res (+ i 1) scale)
(vector-set! res 0 (- (vector-ref res 0)
(* scale (interval-lower-bound domain i))))
(lp (- i 1) coeff)))))))
(define (default-indexer domain)
(coeffs->indexer (default-coeffs domain) domain))
;; Converts the raw integer index to the multi-index in domain that
;; would map to it using the default indexer (i.e. iterating over the
;; possible multi-indices in domain in lexicographic order would
;; produce 0 through volume-1).
(define (invert-default-index domain raw-index)
(let lp ((index raw-index)
(i 0)
(scale (/ (interval-volume domain)
(max 1
(- (interval-upper-bound domain 0)
(interval-lower-bound domain 0)))))
(res '()))
(cond
((>= (+ i 1) (interval-dimension domain))
(reverse (cons (+ index (interval-lower-bound domain i)) res)))
(else
(let ((digit (quotient index scale)))
(lp (- index (* digit scale))
(+ i 1)
(/ scale
(max 1
(- (interval-upper-bound domain (+ i 1))
(interval-lower-bound domain (+ i 1)))))
(cons (+ digit
(interval-lower-bound domain i))
res)))))))
;; Specialized arrays
(define (%make-specialized domain storage body coeffs indexer
safe? mutable? adjacent?)
(%make-array
domain
(specialized-getter body indexer (storage-class-getter storage))
(and mutable?
(specialized-setter body indexer (storage-class-setter storage)))
storage
body
coeffs
indexer
safe?
adjacent?))
(define (make-specialized-array domain . o)
(let* ((storage (if (pair? o) (car o) generic-storage-class))
(safe? (if (and (pair? o) (pair? (cdr o)))
(cadr o)
(specialized-array-default-safe?)))
(body ((storage-class-maker storage)
(interval-volume domain)
(storage-class-default storage)))
(coeffs (default-coeffs domain))
(indexer (coeffs->indexer coeffs domain)))
(assert (boolean? safe?))
(%make-specialized domain storage body coeffs indexer safe? #t #t)))
(define (make-specialized-array-from-data data . o)
(let* ((storage (if (pair? o) (car o) generic-storage-class))
(safe? (if (and (pair? o) (pair? (cdr o)))
(cadr o)
(specialized-array-default-safe?)))
(domain (make-interval (vector ((storage-class-length storage) data))))
(body ((storage-class-data->body storage) data))
(coeffs (default-coeffs domain))
(indexer (coeffs->indexer coeffs domain)))
(assert (boolean? safe?))
(%make-specialized domain storage body coeffs indexer safe? #t #t)))
(define (specialized-array? x)
(and (array? x) (array-storage-class x) #t))
(define (compute-array-packed? array)
(let ((indexer (array-indexer array)))
(call-with-current-continuation
(lambda (return)
(interval-fold
(lambda (prev . multi-index)
(let ((i (apply indexer multi-index)))
(if (and prev (not (= i (+ prev 1))))
(return #f)
i)))
#f
(array-domain array))
#t))))
(define (array-packed? array)
(assert (specialized-array? array))
(let ((res (array-adjacent? array)))
(when (eq? res 'unknown)
(set! res (compute-array-packed? array))
(array-adjacent?-set! array res))
res))
(define (specialized-array-share array new-domain project)
(assert (and (specialized-array? array) (interval? new-domain)))
(let* ((body (array-body array))
(coeffs
(indexer->coeffs
(lambda multi-index
(call-with-values
(lambda () (apply project multi-index))
(array-indexer array)))
new-domain))
(indexer
(coeffs->indexer coeffs new-domain))
(storage (array-storage-class array)))
(%make-specialized new-domain storage body coeffs indexer
(array-safe? array) (array-setter array) 'unknown)))

46
lib/srfi/231/base.sld Normal file
View file

@ -0,0 +1,46 @@
;;> The base array definitions of SRFI 231, plus some extra internal
;;> bindings.
(define-library (srfi 231 base)
(import (scheme base)
(scheme list)
(scheme vector)
(chibi assert))
(export
;; Miscellaneous Functions
translation? permutation?
;; Intervals
make-interval interval? interval-dimension interval-lb interval-ub
interval-lower-bound interval-upper-bound interval-lower-bounds->list
interval-upper-bounds->list interval-lower-bounds->vector
interval-upper-bounds->vector interval= interval-volume
interval-subset? interval-contains-multi-index? interval-projections
interval-for-each interval-dilate interval-intersect
interval-translate interval-permute
interval-scale interval-cartesian-product
interval-width interval-widths
interval-empty?
;; Indexing
index-rotate index-first index-last
indexer->coeffs coeffs->indexer default-indexer default-coeffs
invert-default-index interval-cursor interval-cursor-next!
interval-cursor-next interval-cursor-get interval-fold
;; Storage Classes
make-storage-class storage-class? storage-class-getter
storage-class-setter storage-class-checker storage-class-maker
storage-class-copier storage-class-length storage-class-default
generic-storage-class
storage-class-data? storage-class-data->body
;; Arrays
make-array array? array-domain array-getter array-dimension
mutable-array? array-setter specialized-array-default-safe?
specialized-array-default-mutable?
make-specialized-array make-specialized-array-from-data
specialized-array? array-storage-class array-indexer array-body
array-safe? array-coeffs array-adjacent? array-packed?
specialized-array-share array-ref array-set!
%make-specialized %array-setter-set!
specialized-getter specialized-setter
)
(include "base.scm"))

3825
lib/srfi/231/test.sld Normal file
View file

File diff suppressed because it is too large Load diff

668
lib/srfi/231/transforms.scm Normal file
View file

@ -0,0 +1,668 @@
;; Homogeneous storage classes
;; Define a storage class with an optimized -copy!
(define-syntax define-storage-class
(syntax-rules ()
((define-storage-class name ref set elt? make len default)
(define name
(make-storage-class
ref set elt? make
(lambda (to at from start end)
(let ((limit (min end (+ start (- (len to) at)))))
(if (<= at start)
(do ((i at (+ i 1)) (j start (+ j 1)))
((>= j limit))
(set to i (ref from j)))
(do ((i (+ at (- end start 1)) (- i 1)) (j (- limit 1) (- j 1)))
((< j start))
(set to i (ref from j))))))
len default (lambda (data) #t) (lambda (data) data))))))
(define-storage-class s8-storage-class
s8vector-ref s8vector-set! s8? make-s8vector s8vector-length 0)
(define-storage-class s16-storage-class
s16vector-ref s16vector-set! s16? make-s16vector s16vector-length 0)
(define-storage-class s32-storage-class
s32vector-ref s32vector-set! s32? make-s32vector s32vector-length 0)
(define-storage-class s64-storage-class
s64vector-ref s64vector-set! s64? make-s64vector s64vector-length 0)
(define-storage-class u1-storage-class
u1vector-ref u1vector-set! u1? make-u1vector u1vector-length 0)
(define-storage-class u8-storage-class
u8vector-ref u8vector-set! u8? make-u8vector u8vector-length 0)
(define-storage-class u16-storage-class
u16vector-ref u16vector-set! u16? make-u16vector u16vector-length 0)
(define-storage-class u32-storage-class
u32vector-ref u32vector-set! u32? make-u32vector u32vector-length 0)
(define-storage-class u64-storage-class
u64vector-ref u64vector-set! u64? make-u64vector u64vector-length 0)
(define-storage-class f32-storage-class
f32vector-ref f32vector-set! f32? make-f32vector f32vector-length 0)
(define-storage-class f64-storage-class
f64vector-ref f64vector-set! f64? make-f64vector f64vector-length 0)
(define-storage-class c64-storage-class
c64vector-ref c64vector-set! c64? make-c64vector c64vector-length 0)
(define-storage-class c128-storage-class
c128vector-ref c128vector-set! c128? make-c128vector c128vector-length 0)
(define-storage-class char-storage-class
(lambda (vec i) (integer->char (u32vector-ref vec i)))
(lambda (vec i ch) (u32vector-set! vec i (char->integer ch)))
char? make-u32vector u32vector-length 0)
;; TODO: implement
(define f8-storage-class #f)
(define f16-storage-class #f)
;; Array transformations
(define (array-copy array . o)
(assert (array? array))
(let ((specialized? (specialized-array? array))
(domain (array-domain array)))
(let* ((storage (cond ((pair? o) (car o))
(specialized? (array-storage-class array))
(else generic-storage-class)))
(o (if (pair? o) (cdr o) '()))
(mutable? (cond ((pair? o) (car o))
(specialized? (and (array-setter array) #t))
(else (specialized-array-default-mutable?))))
(o (if (pair? o) (cdr o) '()))
(safe? (cond ((pair? o) (car o))
(specialized? (array-safe? array))
(else (specialized-array-default-safe?)))))
(assert
(and (storage-class? storage) (boolean? mutable?) (boolean? safe?)))
(let* ((body ((storage-class-maker storage)
(interval-volume domain)
(storage-class-default storage)))
(coeffs (default-coeffs domain))
(indexer (coeffs->indexer coeffs domain))
(getter (specialized-getter body indexer
(storage-class-getter storage)))
(setter (specialized-setter body indexer
(storage-class-setter storage)))
(res (%make-specialized domain storage body coeffs indexer
safe? #t #t)))
(array-assign! res array)
(unless mutable?
(%array-setter-set! res #f))
res))))
(define (array-curry array inner-dimension)
(call-with-values
(lambda () (interval-projections (array-domain array) inner-dimension))
(lambda (outer-domain inner-domain)
(cond
((specialized-array? array)
(make-array
outer-domain
(lambda outer-index
(specialized-array-share
array
inner-domain
(lambda inner-index
(apply values (append outer-index inner-index)))))))
(else
(make-array
outer-domain
(lambda outer-index
(make-array
inner-domain
(lambda inner-index
(apply array-ref array (append outer-index inner-index)))
(and
(mutable-array? array)
(lambda (val . inner-index)
(apply array-set! array val (append outer-index inner-index))
))))))))))
(define (array-extract array new-domain)
(assert (and (array? array)
(interval? new-domain)
(interval-subset? new-domain (array-domain array))))
(if (specialized-array? array)
(specialized-array-share array new-domain values)
(make-array new-domain (array-getter array) (array-setter array))))
(define (array-tile array sizes)
(assert (and (array? array)
(vector? sizes)
(= (array-dimension array) (vector-length sizes))
(vector-every exact-integer? sizes)
(vector-every <= sizes (interval-ub (array-domain array)))))
(let ((domain (make-interval
(vector-map
(lambda (lo hi s) (exact (ceiling (/ (- hi lo) s))))
(interval-lb (array-domain array))
(interval-ub (array-domain array))
sizes))))
(make-array
domain
(lambda multi-index
(array-extract
array
(make-interval
(vector-map
(lambda (i lo s) (+ lo (* i s)))
(list->vector multi-index)
(interval-lb (array-domain array))
sizes)
(vector-map
(lambda (i lo hi s)
(min hi (+ lo (* (+ i 1) s))))
(list->vector multi-index)
(interval-lb (array-domain array))
(interval-ub (array-domain array))
sizes)))))))
(define (array-translate array translation)
(let ((new-domain (interval-translate (array-domain array) translation))
(translation-ls (vector->list translation)))
(if (specialized-array? array)
(specialized-array-share
array
new-domain
(lambda multi-index
(apply values (map - multi-index translation-ls))))
(make-array
new-domain
(lambda multi-index
(apply array-ref array (map - multi-index translation-ls)))
(and (mutable-array? array)
(lambda (val . multi-index)
(apply array-set! array val
(map - multi-index translation-ls))))))))
(define (permute ls permutation)
(let ((vec (list->vector ls))
(len (vector-length permutation)))
(do ((i (- len 1) (- i 1))
(res '() (cons (vector-ref vec (vector-ref permutation i)) res)))
((< i 0) res))))
(define (inverse-permutation permutation)
(list->vector
(map car
(list-sort (lambda (a b) (< (cdr a) (cdr b)))
(map cons
(iota (vector-length permutation))
(vector->list permutation))))))
(define (array-permute array permutation)
(assert (permutation? permutation))
(let ((new-domain (interval-permute (array-domain array) permutation))
(perm^-1 (inverse-permutation permutation)))
(if (specialized-array? array)
(specialized-array-share
array
new-domain
(lambda multi-index
(let ((perm-index (permute multi-index perm^-1)))
(apply values perm-index))))
(make-array
new-domain
(lambda multi-index
(let ((perm-index (permute multi-index perm^-1)))
(apply array-ref array perm-index)))
(and (mutable-array? array)
(lambda (val . multi-index)
(apply array-set! array val (permute multi-index perm^-1))))))))
(define (array-reverse array . o)
(assert (array? array))
(let ((flip? (if (pair? o) (car o) (make-vector (array-dimension array) #t))))
(assert (and (vector? flip?)
(= (array-dimension array) (vector-length flip?))
(vector-every boolean? flip?)))
(let* ((flips (vector->list flip?))
(domain (array-domain array))
(lowers (interval-lower-bounds->list domain))
(uppers (interval-upper-bounds->list domain))
(flip-multi-index
(lambda (multi-index)
(map (lambda (i flip-i? lo hi)
(if flip-i? (- (+ lo hi -1) i) i))
multi-index
flips
lowers
uppers))))
(if (specialized-array? array)
(specialized-array-share array
domain
(lambda multi-index
(apply values
(flip-multi-index multi-index))))
(make-array
domain
(lambda multi-index
(apply array-ref array (flip-multi-index multi-index)))
(and
(mutable-array? array)
(lambda (val . multi-index)
(apply array-set! array val (flip-multi-index multi-index))
)))))))
(define (array-sample array scales)
(unless (vector-every zero?
(interval-lower-bounds->vector (array-domain array)))
(error "can only sample an array with zero lower bounds" array))
(let ((scales-ls (vector->list scales))
(new-domain (interval-scale (array-domain array) scales)))
(if (specialized-array? array)
(specialized-array-share
array
new-domain
(lambda multi-index
(apply values (map * multi-index scales-ls))))
(make-array
new-domain
(lambda multi-index
(apply array-ref array (map * multi-index scales-ls)))
(and
(mutable-array? array)
(lambda (val . multi-index)
(apply array-set! array val (map * multi-index scales-ls))))))))
(define (array-outer-product op array1 array2)
(assert (and (procedure? op) (array? array1) (array? array2)))
(make-array (interval-cartesian-product (array-domain array1)
(array-domain array2))
(let ((getter1 (array-getter array1))
(getter2 (array-getter array2))
(dim1 (array-dimension array1)))
(lambda multi-index
(op (apply getter1 (take multi-index dim1))
(apply getter2 (drop multi-index dim1)))))))
(define (array-inner-product A f g B)
(array-outer-product
(lambda (a b) (array-reduce f (array-map g a b)))
(array-copy (array-curry A 1))
(array-copy
(array-curry (array-permute B (index-rotate (array-dimension B) 1))))))
(define (same-dimensions? ls)
(or (null? ls)
(null? (cdr ls))
(and (equal? (array-dimension (car ls)) (array-dimension (cadr ls)))
(same-dimensions? (cdr ls)))))
(define (same-domains? ls)
(or (null? ls)
(null? (cdr ls))
(and (interval= (array-domain (car ls)) (array-domain (cadr ls)))
(same-domains? (cdr ls)))))
(define (array-map f array . arrays)
(make-array (array-domain array)
(let* ((ls (cons array arrays))
(getters (map array-getter ls)))
(assert (same-dimensions? ls))
(lambda multi-index
(apply f (map (lambda (g) (apply g multi-index)) getters))))))
(define (array-for-each f array . arrays)
(if (null? arrays)
(interval-for-each
(let ((g (array-getter array)))
(case (array-dimension array)
((1)
(lambda (i) (f (g i))))
((2)
(lambda (i j) (f (g i j))))
(else
(lambda multi-index
(f (apply g multi-index))))))
(array-domain array))
(interval-for-each
(let* ((lower (interval-lower-bounds->list (array-domain array)))
(ls (cons array arrays))
(getters
(cons (array-getter (car ls))
(map (lambda (ar)
(let ((getter (array-getter ar)))
(lambda multi-index
(apply getter multi-index))))
(cdr ls)))))
(assert (same-domains? ls))
(lambda multi-index
(apply f (map (lambda (g) (apply g multi-index)) getters))))
(array-domain array))))
(define (array-foldl kons knil array)
(interval-fold (lambda (acc . multi-index)
(kons (apply array-ref array multi-index) acc))
knil
(array-domain array)))
(define (array-foldr kons knil array)
(fold-right kons knil (array->list array)))
(define (array-reduce op array)
(let* ((domain (array-domain array))
(init-index (interval-lower-bounds->list domain))
(knil (list 'first-element)))
(interval-fold
(lambda (acc . multi-index)
(if (eq? acc knil)
(apply array-ref array multi-index)
(op acc (apply array-ref array multi-index))))
knil
domain)))
(define (array-any pred array . arrays)
(assert (same-dimensions? (cons array arrays)))
(call-with-current-continuation
(lambda (return)
(apply array-for-each
(lambda args (cond ((apply pred args) => return)))
#f
array
arrays)
#f)))
(define (array-every pred array . arrays)
(assert (same-dimensions? (cons array arrays)))
(call-with-current-continuation
(lambda (return)
(interval-fold
(let ((getters (map array-getter (cons array arrays))))
(lambda (acc . multi-index)
(or (apply pred (map (lambda (g) (apply g multi-index)) getters))
(return #f))))
#t
(array-domain array)))))
(define (array->list array)
(reverse (array-foldl cons '() array)))
(define (list->array domain ls . o)
(let* ((storage (if (pair? o) (car o) generic-storage-class))
(mutable? (if (and (pair? o) (pair? (cdr o)))
(cadr o)
(specialized-array-default-mutable?)))
(safe? (if (and (pair? o) (pair? (cdr o)) (pair? (cddr o)))
(car (cddr o))
(specialized-array-default-safe?)))
(res (make-specialized-array domain storage safe?)))
(assert (and (interval? domain) (storage-class? storage)
(boolean? mutable?) (boolean? safe?)))
(interval-fold
(lambda (ls . multi-index)
(apply array-set! res (car ls) multi-index)
(cdr ls))
ls
domain)
res))
(define (array->vector array)
(list->vector (array->list array)))
(define (vector->array domain vec . o)
(apply list->array domain (vector->list vec o)))
(define (array-assign! destination source)
(assert (and (mutable-array? destination) (array? source)
(interval= (array-domain destination) (array-domain source))))
(let ((getter (array-getter source))
(setter (array-setter destination)))
(interval-for-each
(case (array-dimension destination)
((1) (lambda (i) (setter (getter i) i)))
((2) (lambda (i j) (setter (getter i j) i j)))
((3) (lambda (i j k) (setter (getter i j k) i j k)))
(else
(lambda multi-index
(apply setter (apply getter multi-index) multi-index))))
(array-domain source))
destination))
(define (reshape-without-copy array new-domain)
(let* ((domain (array-domain array))
(orig-indexer (array-indexer array))
(tmp-indexer (default-indexer new-domain))
(new-indexer
(lambda multi-index
(apply orig-indexer
(invert-default-index domain
(apply tmp-indexer multi-index)))))
(new-coeffs (indexer->coeffs new-indexer new-domain #t))
(flat-indexer (coeffs->indexer new-coeffs new-domain))
(new-indexer (coeffs->indexer new-coeffs new-domain))
(body (array-body array))
(storage (array-storage-class array))
(res
(%make-specialized new-domain storage body new-coeffs flat-indexer
(array-safe? array) (array-setter array)
(array-adjacent? array))))
(let ((multi-index (interval-lower-bounds->list domain))
(orig-default-indexer (default-indexer domain)))
(let lp ((i 0)
(ls multi-index))
(let ((reshaped-index
(invert-default-index
new-domain
(apply orig-default-indexer multi-index))))
(cond
((not (equal? (apply flat-indexer reshaped-index)
(apply orig-indexer multi-index)))
#f)
((null? ls)
res)
((= (+ 1 (interval-lower-bound domain i))
(interval-upper-bound domain i))
(lp (+ i 1) (cdr ls)))
(else
(set-car! ls (+ 1 (car ls)))
(lp (+ i 1) (cdr ls)))))))))
(define (specialized-array-reshape array new-domain . o)
(assert (and (specialized-array? array)
(= (interval-volume (array-domain array))
(interval-volume new-domain))))
(let ((copy-on-failure? (and (pair? o) (car o))))
(cond
((reshape-without-copy array new-domain))
(copy-on-failure?
(let ((res (make-specialized-array
new-domain
(array-storage-class array)
(array-safe? array))))
(array-assign! res array)
res))
(else
(error "can't reshape" array new-domain)))))
(define (flatten ls)
(if (pair? (car ls))
(append-map flatten ls)
ls))
(define (list*->array nested-ls . o)
(let lp ((ls nested-ls) (lens '()))
(cond
((pair? ls) (lp (car ls) (cons (length ls) lens)))
(else
(apply list->array
(flatten nested-ls)
(make-interval (list->vector (reverse lens)))
o)))))
(define (array->list* a)
(case (array-dimension a)
((0) (array-ref a))
((1)
(let ((domain (array-domain a)))
(map (lambda (i) (array-ref a i))
(iota (interval-width domain 0)
(interval-lower-bound domain 0)))))
(else
(let ((domain (array-domain a))
(b (array-curry a 1)))
(map (lambda (i) (array->list* (array-ref b i)))
(iota (interval-width domain 0)
(interval-lower-bound domain 0)))))))
(define (vector-iota len start)
(let ((res (make-vector len)))
(do ((i 0 (+ i 1)))
((= i len) res)
(vector-set! res i (+ i start)))))
(define (array->vector* a)
(case (array-dimension a)
((0) (array-ref a))
((1)
(let ((domain (array-domain a)))
(vector-map (lambda (i) (array-ref a i))
(vector-iota (interval-width domain 0)
(interval-lower-bound domain 0)))))
(else
(let ((domain (array-domain a))
(b (array-curry a 1)))
(vector-map (lambda (i) (array->vector* (array-ref b i)))
(vector-iota (interval-width domain 0)
(interval-lower-bound domain 0)))))))
(define (flatten-vec vec)
(if (vector? (vector-ref vec 0))
(append-map flatten-vec vec)
(vector->list vec)))
(define (vector*->array nested-vec . o)
(let lp ((vec nested-vec) (lens '()))
(cond
((vector? vec) (lp (vector-ref vec 0) (cons (vector-length vec) lens)))
(else
(apply list->array
(flatten-vec nested-vec)
(make-interval (list->vector (reverse lens)))
o)))))
(define (dimensions-compatible? a-domain b-domain axis)
(and (= (interval-dimension a-domain) (interval-dimension b-domain))
(let lp ((d (- (interval-dimension a-domain) 1)))
(or (negative? d)
(and (or (= d axis)
(= (- (interval-upper-bound a-domain d)
(interval-lower-bound a-domain d))
(- (interval-upper-bound b-domain d)
(interval-lower-bound b-domain d))))
(lp (- d 1)))))))
(define (array-append axis a . o)
(assert (and (exact-integer? axis)
(array? a)
(< -1 axis (array-dimension a))
(every array? o)))
(let ((a-domain (array-domain a)))
(assert (every (lambda (b)
(dimensions-compatible? a-domain (array-domain b) axis))
o))
(let* ((a-lo (interval-lower-bounds->vector a-domain))
(c-lo (make-vector (interval-dimension a-domain) 0))
(c-hi (interval-widths a-domain)))
(vector-set! c-hi
axis
(fold (lambda (b sum)
(+ sum (interval-width (array-domain b) axis)))
(vector-ref c-hi axis)
o))
(let* ((c-domain (make-interval c-lo c-hi))
(c (make-specialized-array c-domain
(or (array-storage-class a)
generic-storage-class)))
(b-trans (make-vector (array-dimension a) 0)))
(array-assign!
(array-extract c (make-interval c-lo (interval-widths a-domain)))
(array-translate a (vector-map - a-lo)))
(let lp ((arrays o)
(b-offset (- (interval-upper-bound a-domain axis)
(interval-lower-bound a-domain axis))))
(if (null? arrays)
c
(let* ((b (car arrays))
(b-domain (array-domain b))
(b-offset2 (+ b-offset (interval-width b-domain axis)))
(b-lo (make-vector (interval-dimension b-domain) 0))
(b-hi (interval-widths b-domain)))
(vector-set! b-lo axis b-offset)
(vector-set! b-hi axis b-offset2)
(vector-set! b-trans axis (- b-offset))
(let ((view (array-translate
(array-extract c (make-interval b-lo b-hi))
b-trans)))
(array-assign! view b)
(lp (cdr arrays) b-offset2)))))))))
(define (array-stack axis a . o)
(assert (and (exact-integer? axis)
(array? a)
(< -1 axis (array-dimension a))
(every array? o)
(every (lambda (b) (interval= (array-domain a) (array-domain b))) o)))
(let* ((a-lbs (interval-lower-bounds->list (array-domain a)))
(a-ubs (interval-upper-bounds->list (array-domain a)))
(domain
(make-interval
`#(,@(take a-lbs axis) 0 ,@(drop a-lbs axis))
`#(,@(take a-ubs axis) ,(+ 1 (length o)) ,@(drop a-ubs axis))))
(res (make-specialized-array domain
(or (array-storage-class a)
generic-storage-class)))
(perm `#(,axis ,@(delete axis (iota (+ 1 (array-dimension a))))))
(permed (if (zero? axis) res (array-permute res perm)))
(curried (array-curry permed 1))
(get-view (array-getter curried)))
(let lp ((ls (cons a o)) (i 0))
(cond
((null? ls) res)
(else
(array-assign! (get-view i) (car ls))
(lp (cdr ls) (+ i 1)))))))
(define (array-block a . o)
(let ((storage (if (pair? o) (car o) generic-storage-class))
(mutable? (if (and (pair? o) (pair? (cdr o)))
(cadr o)
(specialized-array-default-mutable?)))
(safe? (if (and (pair? o) (pair? (cdr o)) (pair? (cddr o)))
(car (cddr o))
(specialized-array-default-safe?))))
(assert (and (array? a) (not (interval-empty? (array-domain a)))))
(let* ((a-domain (array-domain a))
(get (array-getter a))
(tile0 (apply get (interval-lower-bounds->list a-domain))))
(assert (array? tile0))
(let* ((domain (make-interval
(vector-append (interval-widths a-domain)
(interval-widths (array-domain tile0)))))
(scales (vector->list (interval-widths a-domain)))
(res (make-specialized-array domain storage mutable? safe?)))
(error "TODO: array-block copy data unimplemented")
res))))
(define (array-decurry a . o)
(let ((storage (if (pair? o) (car o) generic-storage-class))
(mutable? (if (and (pair? o) (pair? (cdr o)))
(cadr o)
(specialized-array-default-mutable?)))
(safe? (if (and (pair? o) (pair? (cdr o)) (pair? (cddr o)))
(car (cddr o))
(specialized-array-default-safe?))))
(error "TODO: array-decurry unimplemented")))