| 1 | ;;;;"array.scm" Arrays for Scheme |
| 2 | ; Copyright (C) 1993 Alan Bawden |
| 3 | ; |
| 4 | ; Permission to copy this software, to redistribute it, and to use it |
| 5 | ; for any purpose is granted, subject to the following restrictions and |
| 6 | ; understandings. |
| 7 | ; |
| 8 | ; 1. Any copy made of this software must include this copyright notice |
| 9 | ; in full. |
| 10 | ; |
| 11 | ; 2. Users of this software agree to make their best efforts (a) to |
| 12 | ; return to me any improvements or extensions that they make, so that |
| 13 | ; these may be included in future releases; and (b) to inform me of |
| 14 | ; noteworthy uses of this software. |
| 15 | ; |
| 16 | ; 3. I have made no warrantee or representation that the operation of |
| 17 | ; this software will be error-free, and I am under no obligation to |
| 18 | ; provide any services, by way of maintenance, update, or otherwise. |
| 19 | ; |
| 20 | ; 4. In conjunction with products arising from the use of this material, |
| 21 | ; there shall be no use of my name in any advertising, promotional, or |
| 22 | ; sales literature without prior written consent in each case. |
| 23 | ; |
| 24 | ; Alan Bawden |
| 25 | ; MIT Room NE43-510 |
| 26 | ; 545 Tech. Sq. |
| 27 | ; Cambridge, MA 02139 |
| 28 | ; Alan@LCS.MIT.EDU |
| 29 | |
| 30 | (require 'record) |
| 31 | |
| 32 | ;(declare (usual-integrations)) |
| 33 | |
| 34 | (define array:rtd |
| 35 | (make-record-type "Array" |
| 36 | '(indexer ; Must be a -linear- function! |
| 37 | shape ; Inclusive bounds: ((lower upper) ...) |
| 38 | vector ; The actual contents |
| 39 | ))) |
| 40 | |
| 41 | (define array:indexer (record-accessor array:rtd 'indexer)) |
| 42 | (define array-shape (record-accessor array:rtd 'shape)) |
| 43 | (define array:vector (record-accessor array:rtd 'vector)) |
| 44 | |
| 45 | (define array? (record-predicate array:rtd)) |
| 46 | |
| 47 | (define (array-rank obj) |
| 48 | (if (array? obj) (length (array-shape obj)) 0)) |
| 49 | |
| 50 | (define (array-dimensions ra) |
| 51 | (map (lambda (ind) (if (zero? (car ind)) (+ 1 (cadr ind)) ind)) |
| 52 | (array-shape ra))) |
| 53 | |
| 54 | (define array:construct |
| 55 | (record-constructor array:rtd '(shape vector indexer))) |
| 56 | |
| 57 | (define (array:compute-shape specs) |
| 58 | (map (lambda (spec) |
| 59 | (cond ((and (integer? spec) |
| 60 | (< 0 spec)) |
| 61 | (list 0 (- spec 1))) |
| 62 | ((and (pair? spec) |
| 63 | (pair? (cdr spec)) |
| 64 | (null? (cddr spec)) |
| 65 | (integer? (car spec)) |
| 66 | (integer? (cadr spec)) |
| 67 | (<= (car spec) (cadr spec))) |
| 68 | spec) |
| 69 | (else (slib:error "array: Bad array dimension: " spec)))) |
| 70 | specs)) |
| 71 | |
| 72 | (define (make-array initial-value . specs) |
| 73 | (let ((shape (array:compute-shape specs))) |
| 74 | (let loop ((size 1) |
| 75 | (indexer (lambda () 0)) |
| 76 | (l (reverse shape))) |
| 77 | (if (null? l) |
| 78 | (array:construct shape |
| 79 | (make-vector size initial-value) |
| 80 | (array:optimize-linear-function indexer shape)) |
| 81 | (loop (* size (+ 1 (- (cadar l) (caar l)))) |
| 82 | (lambda (first-index . rest-of-indices) |
| 83 | (+ (* size (- first-index (caar l))) |
| 84 | (apply indexer rest-of-indices))) |
| 85 | (cdr l)))))) |
| 86 | |
| 87 | (define (make-shared-array array mapping . specs) |
| 88 | (let ((new-shape (array:compute-shape specs)) |
| 89 | (old-indexer (array:indexer array))) |
| 90 | (let check ((indices '()) |
| 91 | (bounds (reverse new-shape))) |
| 92 | (cond ((null? bounds) |
| 93 | (array:check-bounds array (apply mapping indices))) |
| 94 | (else |
| 95 | (check (cons (caar bounds) indices) (cdr bounds)) |
| 96 | (check (cons (cadar bounds) indices) (cdr bounds))))) |
| 97 | (array:construct new-shape |
| 98 | (array:vector array) |
| 99 | (array:optimize-linear-function |
| 100 | (lambda indices |
| 101 | (apply old-indexer (apply mapping indices))) |
| 102 | new-shape)))) |
| 103 | |
| 104 | (define (array:in-bounds? array indices) |
| 105 | (let loop ((indices indices) |
| 106 | (shape (array-shape array))) |
| 107 | (if (null? indices) |
| 108 | (null? shape) |
| 109 | (let ((index (car indices))) |
| 110 | (and (not (null? shape)) |
| 111 | (integer? index) |
| 112 | (<= (caar shape) index (cadar shape)) |
| 113 | (loop (cdr indices) (cdr shape))))))) |
| 114 | |
| 115 | (define (array:check-bounds array indices) |
| 116 | (or (array:in-bounds? array indices) |
| 117 | (slib:error "array: Bad indices for " array indices))) |
| 118 | |
| 119 | (define (array-ref array . indices) |
| 120 | (array:check-bounds array indices) |
| 121 | (vector-ref (array:vector array) |
| 122 | (apply (array:indexer array) indices))) |
| 123 | |
| 124 | (define (array-set! array new-value . indices) |
| 125 | (array:check-bounds array indices) |
| 126 | (vector-set! (array:vector array) |
| 127 | (apply (array:indexer array) indices) |
| 128 | new-value)) |
| 129 | |
| 130 | (define (array-in-bounds? array . indices) |
| 131 | (array:in-bounds? array indices)) |
| 132 | |
| 133 | ; Fast versions of ARRAY-REF and ARRAY-SET! that do no error checking, |
| 134 | ; and don't cons intermediate lists of indices: |
| 135 | |
| 136 | (define (array-1d-ref a i0) |
| 137 | (vector-ref (array:vector a) ((array:indexer a) i0))) |
| 138 | |
| 139 | (define (array-2d-ref a i0 i1) |
| 140 | (vector-ref (array:vector a) ((array:indexer a) i0 i1))) |
| 141 | |
| 142 | (define (array-3d-ref a i0 i1 i2) |
| 143 | (vector-ref (array:vector a) ((array:indexer a) i0 i1 i2))) |
| 144 | |
| 145 | (define (array-1d-set! a v i0) |
| 146 | (vector-set! (array:vector a) ((array:indexer a) i0) v)) |
| 147 | |
| 148 | (define (array-2d-set! a v i0 i1) |
| 149 | (vector-set! (array:vector a) ((array:indexer a) i0 i1) v)) |
| 150 | |
| 151 | (define (array-3d-set! a v i0 i1 i2) |
| 152 | (vector-set! (array:vector a) ((array:indexer a) i0 i1 i2) v)) |
| 153 | |
| 154 | ; STOP! Do not read beyond this point on your first reading of |
| 155 | ; this code -- you should simply assume that the rest of this file |
| 156 | ; contains only the following single definition: |
| 157 | ; |
| 158 | ; (define (array:optimize-linear-function f l) f) |
| 159 | ; |
| 160 | ; Of course everything would be pretty inefficient if this were really the |
| 161 | ; case, but it isn't. The following code takes advantage of the fact that |
| 162 | ; you can learn everything there is to know from a linear function by |
| 163 | ; simply probing around in its domain and observing its values -- then a |
| 164 | ; more efficient equivalent can be constructed. |
| 165 | |
| 166 | (define (array:optimize-linear-function f l) |
| 167 | (let ((d (length l))) |
| 168 | (cond |
| 169 | ((= d 0) |
| 170 | (array:0d-c (f))) |
| 171 | ((= d 1) |
| 172 | (let ((c (f 0))) |
| 173 | (array:1d-c0 c (- (f 1) c)))) |
| 174 | ((= d 2) |
| 175 | (let ((c (f 0 0))) |
| 176 | (array:2d-c01 c (- (f 1 0) c) (- (f 0 1) c)))) |
| 177 | ((= d 3) |
| 178 | (let ((c (f 0 0 0))) |
| 179 | (array:3d-c012 c (- (f 1 0 0) c) (- (f 0 1 0) c) (- (f 0 0 1) c)))) |
| 180 | (else |
| 181 | (let* ((v (map (lambda (x) 0) l)) |
| 182 | (c (apply f v))) |
| 183 | (let loop ((p v) |
| 184 | (old-val c) |
| 185 | (coefs '())) |
| 186 | (cond ((null? p) |
| 187 | (array:Nd-c* c (reverse coefs))) |
| 188 | (else |
| 189 | (set-car! p 1) |
| 190 | (let ((new-val (apply f v))) |
| 191 | (loop (cdr p) |
| 192 | new-val |
| 193 | (cons (- new-val old-val) coefs))))))))))) |
| 194 | |
| 195 | ; 0D cases: |
| 196 | |
| 197 | (define (array:0d-c c) |
| 198 | (lambda () c)) |
| 199 | |
| 200 | ; 1D cases: |
| 201 | |
| 202 | (define (array:1d-c c) |
| 203 | (lambda (i0) (+ c i0))) |
| 204 | |
| 205 | (define (array:1d-0 n0) |
| 206 | (cond ((= 1 n0) +) |
| 207 | (else (lambda (i0) (* n0 i0))))) |
| 208 | |
| 209 | (define (array:1d-c0 c n0) |
| 210 | (cond ((= 0 c) (array:1d-0 n0)) |
| 211 | ((= 1 n0) (array:1d-c c)) |
| 212 | (else (lambda (i0) (+ c (* n0 i0)))))) |
| 213 | |
| 214 | ; 2D cases: |
| 215 | |
| 216 | (define (array:2d-0 n0) |
| 217 | (lambda (i0 i1) (+ (* n0 i0) i1))) |
| 218 | |
| 219 | (define (array:2d-1 n1) |
| 220 | (lambda (i0 i1) (+ i0 (* n1 i1)))) |
| 221 | |
| 222 | (define (array:2d-c0 c n0) |
| 223 | (lambda (i0 i1) (+ c (* n0 i0) i1))) |
| 224 | |
| 225 | (define (array:2d-c1 c n1) |
| 226 | (lambda (i0 i1) (+ c i0 (* n1 i1)))) |
| 227 | |
| 228 | (define (array:2d-01 n0 n1) |
| 229 | (cond ((= 1 n0) (array:2d-1 n1)) |
| 230 | ((= 1 n1) (array:2d-0 n0)) |
| 231 | (else (lambda (i0 i1) (+ (* n0 i0) (* n1 i1)))))) |
| 232 | |
| 233 | (define (array:2d-c01 c n0 n1) |
| 234 | (cond ((= 0 c) (array:2d-01 n0 n1)) |
| 235 | ((= 1 n0) (array:2d-c1 c n1)) |
| 236 | ((= 1 n1) (array:2d-c0 c n0)) |
| 237 | (else (lambda (i0 i1) (+ c (* n0 i0) (* n1 i1)))))) |
| 238 | |
| 239 | ; 3D cases: |
| 240 | |
| 241 | (define (array:3d-01 n0 n1) |
| 242 | (lambda (i0 i1 i2) (+ (* n0 i0) (* n1 i1) i2))) |
| 243 | |
| 244 | (define (array:3d-02 n0 n2) |
| 245 | (lambda (i0 i1 i2) (+ (* n0 i0) i1 (* n2 i2)))) |
| 246 | |
| 247 | (define (array:3d-12 n1 n2) |
| 248 | (lambda (i0 i1 i2) (+ i0 (* n1 i1) (* n2 i2)))) |
| 249 | |
| 250 | (define (array:3d-c12 c n1 n2) |
| 251 | (lambda (i0 i1 i2) (+ c i0 (* n1 i1) (* n2 i2)))) |
| 252 | |
| 253 | (define (array:3d-c02 c n0 n2) |
| 254 | (lambda (i0 i1 i2) (+ c (* n0 i0) i1 (* n2 i2)))) |
| 255 | |
| 256 | (define (array:3d-c01 c n0 n1) |
| 257 | (lambda (i0 i1 i2) (+ c (* n0 i0) (* n1 i1) i2))) |
| 258 | |
| 259 | (define (array:3d-012 n0 n1 n2) |
| 260 | (cond ((= 1 n0) (array:3d-12 n1 n2)) |
| 261 | ((= 1 n1) (array:3d-02 n0 n2)) |
| 262 | ((= 1 n2) (array:3d-01 n0 n1)) |
| 263 | (else (lambda (i0 i1 i2) (+ (* n0 i0) (* n1 i1) (* n2 i2)))))) |
| 264 | |
| 265 | (define (array:3d-c012 c n0 n1 n2) |
| 266 | (cond ((= 0 c) (array:3d-012 n0 n1 n2)) |
| 267 | ((= 1 n0) (array:3d-c12 c n1 n2)) |
| 268 | ((= 1 n1) (array:3d-c02 c n0 n2)) |
| 269 | ((= 1 n2) (array:3d-c01 c n0 n1)) |
| 270 | (else (lambda (i0 i1 i2) (+ c (* n0 i0) (* n1 i1) (* n2 i2)))))) |
| 271 | |
| 272 | ; ND cases: |
| 273 | |
| 274 | (define (array:Nd-* coefs) |
| 275 | (lambda indices (apply + (map * coefs indices)))) |
| 276 | |
| 277 | (define (array:Nd-c* c coefs) |
| 278 | (cond ((= 0 c) (array:Nd-* coefs)) |
| 279 | (else (lambda indices (apply + c (map * coefs indices)))))) |