[clean-list] Recursive array algorithms

[Mike Rainey]

Hi, I'm David Wise's student, and am attempting to implement some recursive 
array libraries in Clean.  In particular, I'm adding support for quadtree 
matrices.  The timing results will (hopefully) generate comparison graphs 
against Haskell and C, which will later end up in a journal submission.  
Clean, with its uniqueness type system, seems an excellent environment for 
quadtree matrices.  However, there are some fundamental difficulties that 
preclude good performance.  To explain my problem I'll start with the 
canonical example, matrix multiplication.

matrMult::*{#Real} {#Real} {#Real} !Int -> *{#Real} 
matrMult c a b ord = mm c (MkM 0) (MkM 0) (MkM 0) 0
	where 
	dpth = lg2 ord   // matrix depth (in number of recursive steps to reach the 
base case) for an n x n size is lg base 2 of n
	mm::*{#Real} !MortonIx !MortonIx !MortonIx !Int -> *{#Real}
	mm c cin ain bin lvl | lvl >= dpth =  // base case (1x1 block)
							let (cv,cp) = c![coerceInt(cin)] in  
								{cp & [coerceInt(cin)]= cv + a.[coerceInt(ain)] * b.[coerceInt(bin)]}
   			              | otherwise #     c = mm c (nw cin) (nw ain) (nw bin) lvlp
	   				        	       c = mm c (nw cin) (ne ain) (sw bin) lvlp
							       c = mm c (ne cin) (nw ain) (ne bin) lvlp
							       c = mm c (ne cin) (ne ain) (se bin) lvlp
						               c = mm c (sw cin) (sw ain) (nw bin) lvlp
							       c = mm c (sw cin) (se ain) (sw bin) lvlp
							       c = mm c (se cin) (sw ain) (ne bin) lvlp
					                 =         mm c (se cin) (se ain) (se bin) lvlp
	   where lvlp = lvl + 1

matrMult takes matrices c, a, and b and stores the result of a*b into c.  The 
mm helper function multiplies the 1x1 base case in the first guard, and 
recursively computes the n x n case in the second guard.  To recurse down 
this matrix, Morton order indices are used.  They are selected by quadtant as 
nw (northwest), sw (south west), etc.  This test code has no bounds checking, 
as it is only intended for powers of 2.

The point is, however, that the mm function, with it's 8 recursive calls, ends 
up clobbering the stack.  The program segfaults for interesting sizes (those 
> 64x64).  Upping the stack size solves this, until size is greater than 
2048x2048.  But this is a poor solution, as a constant size stack will do.  
Should Clean recognize that mm is tail-recursive?  Is there some magic 
incantation that might force mm to be tail-recursive?  I'm very interested in 
how uniqueness types can increase performance, and would appreciate any 
insights.

Thanks!
Mike Rainey
Indiana University CS, soon to be University of Chicago CS