<div>On 8 April 2012 18:59, Caligo <span dir="ltr"><<a href="mailto:iteronvexor@gmail.com">iteronvexor@gmail.com</a>></span> wrote:</div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div class="im">On Tue, Apr 3, 2012 at 6:20 AM, Cristi Cobzarenco<br>
<<a href="mailto:cristi.cobzarenco@gmail.com">cristi.cobzarenco@gmail.com</a>> wrote:<br>
><br>
> The point of these is to have light-weight element wise operation support.<br>
> It's true that in theory the built-in arrays do this. However, this library<br>
> is built on top BLAS/LAPACK, which means operations on large matrices will<br>
> be faster than on D arrays.<br>
><br>
<br>
</div>I can't agree with building it on top of LAPACK or any other BLAS<br>
implementation, but perhaps I shouldn't complain because I'm not the<br>
one who's going to be implementing it. I like the approach Eigen has<br>
taken where it offers its own BLAS implementation and, iirc, other<br>
BLAS libraries can be used as optional back-ends.<br></blockquote><div><br></div><div>Yes, I agree with you. I already built naive implementations for BLAS & LAPACK functions as part of my last year project, using external libraries is optional (building with version=nodeps ensures absolutely no dependencies are needed). The argument still stands, though. If you _do_ use external BLAS libraries then using value arrays will _still_ be faster.</div>
<div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div class="im"><br>
> Also, as far as I know, D doesn't support<br>
> allocating dynamic 2-D arrays (as in not arrays of arrays), not to mention<br>
> 2-D slicing (keeping track of leading dimension).<br>
><br>
<br>
</div>I fail to see why there is any need for 2D arrays. We need to make<br>
sure multidimensional arrays (matrices) have data in very good<br>
arrangements. This is called tiling and it requires 1D arrays: 2D<br>
arrays are stored as 1D arrays together with an indexing mechanism.<br></blockquote><div><br></div><div>That is precisely what I meant. We need wrappers around D arrays, because they, by themselves, do not support 2-D indexing. By providing a wrapper we allow the nice syntax of matrix[ a, b ]. This kind of wrapping is already in SciD by CowMatrix. I meant we shouldn't use D built-in arrays directly, not not at all. Also, as mentioned before, we can't use new for allocation, because we want the library to be GC-independent.</div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div class="im"><br>
> Also I'm not sure how a case like this will be compiled, it may or may not<br>
> allocate a temporary:<br>
><br>
> a[] = b[] * c[] + d[] * 2.0;<br>
><br>
> The expression templates in SciD mean there will be no temporary allocation<br>
> in this call.<br>
><br>
<br>
</div>Why are expression templates used? <br></blockquote><div><br></div><div>As H. S. Teoh rightfully pointed out, it is important not to allocate temporaries in matrix operations. You want to evaluate A = B * 3.0 + D * 2.0 (where .* is element-wise multiplication) as (in BLAS terms):</div>
<div> copy( B, A );</div><div> scal( 3.0, A );</div><div> axpy( D, 2.0, A );</div><div><br></div><div>Or with naive expression template evaluation:</div><div> for( i = 0 ; i < n ; ++ i ) {</div><div> A[ i ] = B[ i ] * 3.0 + D * 2.0;</div>
<div> }</div><div><br></div><div>The D compiler would instead evaluate this as:</div><div> // allocate temporaries</div><div> allocate( tmp1, A.length );</div><div> allocate( tmp2, A.length );</div><div> allocate( tmp3, A.length );</div>
<div><br></div><div> // compute tmp1 = B * 3.0</div><div> copy( B, tmp1 );</div><div> scal( 3.0, tmp1 );</div><div><br></div><div> // compute tmp2 = D * 2.0;</div><div> copy( D, tmp2 );</div><div> scal( 2.0, tmp2 );</div>
<div><br></div><div> // compute tmp3 = tmp1 + tmp2;</div><div> copy( tmp1, tmp3 );</div><div> axpy( tmp2, 1.0, tmp1 );</div><div><br></div><div> // copy tmp3 into A</div><div> copy( tmp3, A );</div><div><br></div><div>
Plenty of useless computation. Note this is not a fault of the compiler, it has no way of knowing which temporaries can be optimised away for user types (i.e. our matrices).</div><div><br></div><div>> Are you pretty much rewriting Eigen in D?</div>
<div><br></div><div>No. It is just the interface - and even that only to a certain extent - that will mimic Eigen. The core the library would be very similar to what I implemented for SciD last year, which, you will find, is very D-like and not at all like Eigen.</div>
<div><br></div><div><br></div></div>