SRFI 231: array-inner-product fix

[gambiteer]

You have

> chibi-scheme 
> (import (srfi 231))
> (array-inner-product (make-array (make-interval '#(4 0)) list)
                     list
                     list
                     (make-array (make-interval '#(0 4)) list))
ERROR in array-inner-product on line 332 of file /usr/local/chibi//share/chibi/srfi/231/transforms.scm: not enough args
    #<procedure array-curry 2>
    1

Fixed with:

 git diff
diff --git a/lib/srfi/231/transforms.scm b/lib/srfi/231/transforms.scm
index e3a97955..2d6802f9 100644
--- a/lib/srfi/231/transforms.scm
+++ b/lib/srfi/231/transforms.scm
@@ -325,11 +325,21 @@
                       (apply getter2 (drop multi-index dim1)))))))

 (define (array-inner-product A f g B)
+  (assert (and (array? A) (array? B)
+               (procedure? f) (procedure? g)
+               (positive? (array-dimension A)) (positive? (array-dimension B))
+               (let ((A-dim (array-dimension A))
+                     (A-dom (array-domain A))
+                     (B-dom (array-domain B)))
+                 (and (eqv? (interval-lower-bound A-dom (- A-dim 1))
+                            (interval-lower-bound B-dom 0))
+                      (eqv? (interval-upper-bound A-dom (- A-dim 1))
+                            (interval-upper-bound B-dom 0))))))
   (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))))))
+    (array-curry (array-permute B (index-rotate (array-dimension B) 1)) 1))))

 (define (same-dimensions? ls)
   (or (null? ls)

(along with some error checking).

[ashinn]

Thanks, but your fix results in a lazily evaluated array-reduce on an empty array, so for instance calling array->list* on the result gives the error:

    ERROR: (array->list* (array-inner-product (make-array (make-interval '...
        Exception: assertion failed
    (not (interval-empty? domain))
    (domain #<Interval 134075314546304>)
        on line 3589 of file "./lib/srfi/231/test.sld"
        (array->list* (array-inner-product (make-array (make-interval '#(4 0)) list) list list (make-array (make-interval '#(0 4)) list)))

What is the correct semantics here?

[gambiteer]

That's a good question. You're making something from nothing, basically.

It appears that APL does not allow empty arrays (intrinsically).

NumPy seems to suffer from the same problem: https://stackoverflow.com/questions/45282176/numpy-dot-product-with-inner-dimension-zero-gives-unexpected-result

In Matlab, matrix products are the only type of inner product, it appears, so there are multiplicative and additive identities, so reduce gives a reasonable answer: https://www.mathworks.com/help/matlab/math/empty-arrays.html See also: https://www.mathworks.com/help/matlab/matlab_prog/compatible-array-sizes-for-basic-operations.html

So I would say "It is an error if the width of the last axis of A (the same as the first axis of B) is zero".

Then it would just be a question of where to catch that error.

I'd be inclined to add (not (zero? (interval-width B-dom 0))) to the second (and ...) in the assert.

[gambiteer]

It seems that the sample implementation recognized that this situation could happen; accessing an entry of the result of array-inner-product threw an error with a message beginning array-inner-product: Attempting to reduce over an empty array: I.e., the name of the caller (array-inner-product) was passed to the internal procedure %%array-reduce for the error message.

So the sample implementation delayed the error until there was an attempt to access an entry of the result of array-inner-product.

I'm not sure I agree with that point of view now.

[mnieper]

In Matlab, matrix products are the only type of inner product, it appears, so there are multiplicative and additive identities, so reduce gives a reasonable answer: https://www.mathworks.com/help/matlab/math/empty-arrays.html See also: https://www.mathworks.com/help/matlab/matlab_prog/compatible-array-sizes-for-basic-operations.html

So I would say "It is an error if the width of the last axis of A (the same as the first axis of B) is zero".

Unless you have a neutral element for the addition (which you don't supply array-inner-product with), it has to be an error.

In most natural applications, however, there is a neutral element, so array-inner-product does not catch the common situations quite right. A better version would take a ring (e.g., modelled by a record) that encodes addition, multiplication, and neutral elements.

Then it would just be a question of where to catch that error.

I'd be inclined to add (not (zero? (interval-width B-dom 0))) to the second (and ...) in the assert.

Dimensions are never lazily calculated, are they? If so, the check should be performed eagerly.

[ashinn]

Rather than introducing a ring concept, taking a neutral element as an optional parameter is easier, but the SRFI is finalized. Changed to fail fast.

[gambiteer]

A better version would take a ring (e.g., modelled by a record) that encodes addition, multiplication, and neutral elements.

I want to push back on this a bit.

The data in arrays have no intrinsic algebraic or mathematical structure. I haven't wanted to complicate this SRFI (it's complicated enough!) by building in the possibility of specifying algebraic structures for some of the operations.