[Haskell] Eliminating Multiple-Array Bound Checking through Non-dependent types

An earlier message showed an example of writing code with non-trivial static guarantees in the present-day Haskell (i.e., Haskell98 + rank-2 types). http://pobox.com/~oleg/ftp/Haskell/types.html#branding Although this approach obviously falls short of the genuine dependent typing, we can use it in the existing, well-maintained functional language systems. Thus we may begin to acquire taste for dependent typing -- which, McBride at al. point out, matter. The present message shows a more advanced example -- eliminating array bound checking when processing several arrays at a time. The number of arrays to process is not statically known. Furthermore, the arrays may have different sizes and bounds -- potentially, empty and non-overlapping bounds. The example involves native Haskell arrays, index computations, and general recursion. All array indexing operations are statically guaranteed to be safe -- and so we can safely use an efficient unsafeAt provided by GHC seemingly for that purpose. The code is efficient; the static assurances in the main loop over arrays cost us no run-time overhead. The example uses only Haskell98 + higher-ranked types. No new type classes are introduced. The safety is based on: Haskell type system, quantified type variables, and a compact general-purpose trusted kernel. I thank Daniel Yokomizo for the challenge. Our example is folding over multiple, variously-sized arrays. This is like a fold over an array -- generalized to an arbitrary number of arrays, whose index ranges do not have to be the same (and do not have to overlap). Typing this example in a genuinely dependent type system is probably going to be quite challenging. Our goal is to implement a function marray_fold :: (Ix i, Integral i) => (seed -> [e] -> seed) -> seed -> [Array i e] -> seed Its third argument is a list of arrays; the arrays have all the same element and index types; the actual sizes (that is, the lower and upper bounds) may differ. Some arrays in the list may even be empty (with the lower bound higher than the upper one). The function marray_fold, like left fold, applies its left argument to the values extracted from the corresponding arrays. Because arrays may have different sizes and bounds, marray_fold operates over the range of indices that is the intersection of the ranges of the argument arrays. For example: dot a1 a2 = marray_fold (\seed l -> product l + seed) 0 [a1,a2] computes the dot products of two arrays. Granted our example is a bit contrived: in practice one would probably use ``fixed arity'' multiple array fold (similar to zipWith, zipWith3, etc -- which are more efficient and simpler to implement). Also simple is the processing of arrays of the same bounds (and invoking an error continuations if the bounds turn out different). OTH, our present example is more challenging and thus more interesting. Here's the complete implementation (assuming the short trusted branding library explained below). > {-# OPTIONS -fglasgow-exts #-} > module MDep where > import Data.Array > > marray_fold :: (Ix i, Integral i) => > (seed -> [e] -> seed) -> seed -> [Array i e] -> seed > marray_fold _ seed [] = seed -- No arrays, no fold > marray_fold f seed arrs = brands arrs (marray_fold' f seed) seed > > marray_fold' f seed (barrs, bl, bh) = loop (bl2i bl) seed > where loop bi seed = bsucc bh bi (\bi' -> loop bi' seed') seed' > where seed' = f seed (map (!. bi) barrs) > > brands :: (Ix i, Integral i) => [Array i e] > -> (forall s. ([BArray s i e], BLow s i, BHi s i) -> w) > -> w -> w > brands [arr] consumer onempty = > brand arr (\ (barr,bl,bh) -> consumer ([barr],bl,bh)) onempty > brands (a:arrs) consumer onempty = > brands arrs (\bbars -> brand_merge bbars a consumer onempty) onempty The function |brands| computes the intersection index range for the given arrays, and brands each array and the intersection range. The range is given as a pair (BLow s i, BHi s i). The quantified (and thus unforgeable) type variable |s| represents the brand. The function |brands| has a higher-ranked type and must have a signature. Incidentally, functions like that are relatively rare; processors of the list of arrays of the same brand do not require any signature (cf. marray_fold'). The folding loop has as few index comparison operations as algorithmically needed, i.e., only one per iteration. Operator (!.) is a statically safe array indexing operator. The type system and the trust properties of the kernel below guarantee that in the expression "arr !. m" the index `m' is positively in range of the array `arr' bounds. Here are a few test cases, of the dot product over several arrays > dots ars = marray_fold (\seed l -> product l + seed) 0 ars > test1 = dots [listArray (1,5) [1..]] > test2 = dots [listArray (1,5) [1..],listArray (1,5) [1,1..]] > test3 = dots [listArray (1,5) [1..],listArray (5,1) [1,1..]] > test4 = dots [listArray (1,5) [1..],listArray (2,7) [1,1..]] > test5 = dots [listArray (1,5) [1..],listArray (7,7) [1,1..]] > test6 = dots [listArray (1,5) [1..],listArray (5,5) [1,1..], > listArray (3,7) [3..]] The case test3 demonstrates the inclusion of an empty array; test5 shows non-overlapping arrays. The code relies on the compact general-purpose trusted kernel explained below. That code should be put into a separate module. First we introduce tags for Branded arrays, Branded indices, and Branded low and high bounds. > -- the data constructors of those four must *not* be exported! > newtype BArray s i a = BArray (Array i a) > newtype BIndex s i = BIndex i > newtype BLow s i = BLow i > newtype BHi s i = BHi i > unbi (BIndex i) = i These are `newtype's and so impose no run-time overhead. Of interest is a phantom type variable 's', which marks a _brand_ of an array and of an array index. An index or a bound are branded if they are positively within the range of the array of their brand. The type variable 's' is similar to that in the ST monad. The latter relies on 's' to enforce serialization. We, OTH, do not impose any linearity constraints on 's' -- it may be freely duplicated and discarded (see `unbi'). It is created however under controlled conditions. The safety depends on the trusted way of creating branded types: the data constructors BIndex, BLow, BHi, and BArray must be used in the trusted kernel only, and should not be available anywhere else. The uniqueness of 's' afforded by the explicit universal quantification prevents mixing up of different brands. We must re-iterate that safety depends on assurances of the code that constructs branded values. Because of the high assurance, we must formulate the safety properties as propositions, and prove them. Fortunately, the code below is compact and straightforward, as well as general purpose. We start with an introduction rule for BArray: > brand:: (Ix i, Integral i) > => Array i e > -> (forall s. (BArray s i e, BLow s i, BHi s i) -> w) > -> w -> w > brand (a::Array i e) k kempty = > let (l,h) = bounds a > in if l <= h then k ((BArray a)::BArray () i e, BLow l, BHi h) > else kempty The function takes an array and two continuations. If the array is empty, kempty is returned. Otherwise, the continuation k is invoked, which receives the branded array along with two branded indices, for the lower and the upper bounds of the array. Proposition: For the argument (a,l,h) of the continuation k, the branded integral indices l and h are within the bounds of the non-empty branded array a. Proof: from the semantics of the functions `bounds' and `inRange', properties of integral numbers and ranges over integral numbers (specifically, contiguity), and the code itself. The function brand has a higher-rank type. It is the existential quantification of 's' as well as the absence of BArray constructor elsewhere guarantee that the same brand entails the same bounds. For this example, we need another introduction rule for BArray: > brand_merge :: (Ix i, Integral i) => > ([BArray s i e], BLow s i, BHi s i) > -> Array i e > -> (forall s'. ([BArray s' i e], BLow s' i, BHi s' i) -> w) > -> w -> w > brand_merge (bas,(BLow bl),(BHi bh)) (a :: Array i e) k kempty = > let (l,h) = bounds a > l' = max l bl > h' = min h bh > in if l' <= h' then > k (((BArray a)::BArray () i e) : > (map (\ (BArray a) -> BArray a) bas), > BLow l', BHi h') > else kempty The function |brand_merge| takes a list of branded arrays and their greatest common bound (specified as a pair of branded indices). The function also takes an unbranded array. We prepend the new array to the existing ones, and adjust the greatest common bound. The result will have a new brand, which signifies that the new |bl| and |bh| indices are within the range of all of the arrays. The mapped function (\ (BArray a) -> BArray a) is operationally identity. Hopefully the compiler recognizes this fact. Proposition: Given (bas, l, h) where [l,h] is an integral range, forall arr in bas. [l,h] in bounds(arr) if the continuation k is ever invoked, its argument (bas',l',h') has the property: bas' is a non-empty list, [l',h'] is a non-empty linear (integral) range, forall arr in bas'. [l',h'] in bounds(arr) Proof: from the semantics of the functions `bounds' and `inRange', properties of integral numbers and ranges over integral numbers (specifically, contiguity), and the code itself. > bl2i:: BLow s i -> BIndex s i > bl2i (BLow i) = BIndex i Proposition: if bl is within the bounds of an array branded with s, then the result is BIndex within the bounds of that array. Proof: trivial. > bsucc:: (Integral i, Ix i) > => BHi s i -> BIndex s i -> (BIndex s i -> r) -> r -> r > bsucc (BHi upb) (BIndex i) on_within on_out > = let i' = succ i > in if i' <= upb then (on_within (BIndex i')) else on_out The function `bsucc' takes two branded indices that correspond to the same bounds (see the variable 's'). The function also takes two continuations, on_within and on_out. The first index is considered to be an upper limit. The function increments the second index. If the result does not exceed the upper limit, we invoke on_within and pass it the result. Otherwise, we invoke `on_out'. Safety Proposition: l <= upb <= h, l <= i <= h, (i+1) <= upb |- l <= (i+1) <= h Proof: from i < (i+1) and properties of inequalities. The safety proposition justifies our use of the data constructor BIndex. > bpred:: (Ix i, Integral i) > => BLow s i -> BIndex s i -> (BIndex s i -> r) -> r -> r > bpred (BLow lwb) (BIndex i) on_within on_out > = let i' = pred i > in if i' >= lwb then (on_within (BIndex i')) else on_out The dual of `bsucc'. Safety Proposition: l <= lwb <= h, l <= i <= h, (i-1) >= lwb |- l <= (i-1) <= h Because a branded index is assuredly within the bounds of the array of the same brand, we can write > infixl 5 !. > (!.):: (Ix i) => BArray s i e -> (BIndex s i) -> e > (BArray a) !. (BIndex i) = a ! i actually, we may _safely_ replace a ! i with `unsafeAt a i' A reader may have noticed a seemingly artificial |Integral i| constraint on indices. Why can't the mere |Ix i| constraint suffice? That would let us potentially deal with multi-dimensional arrays. Alas, the sole |Ix i| constraint does not suffice because the Ix laws guaranteed in Haskell are lame. For example: inRange (l,h) l ? Does not always hold: inRange (2,1) 2 ---> False if rangeSize (l,h) >0, inRange (l,h) l ? does not seem to follow from the provided laws l <= h, inRange(l,h) l ? False: let l = (1,2); h = (2,1) in (l<=h, inRange (l,h) l) gives (True,False) let (bl,bh) = ((1,2),(1,2)) -- non-empty range (l,h) = ((1,1),(2,1)) -- non-empty range l' = max l bl -- (1,2) h' = min h bh -- (1,2) rangeSize (l',h') is positive but inRange (l,h) l' is False Although Ord is a super-class of Ix, the comparison operation of Ord is not required to have any relationship with ranges. And in fact, it doesn't: the comparison operation on tuples is lexicographic, which isn't the best choice for tuples representing multi-dimensional ranges. Perhaps this message might provide the incentive for more precise specification of the laws for various basic classes in future Haskell standards.