ch-ch-update: series, closed-form series, and strides

[Bill Baxter]

On Sat, Jan 31, 2009 at 1:40 PM, Andrei Alexandrescu
<SeeWebsiteForEmail at erdani.org> wrote:
> Bill Baxter wrote:
>>
>> On Sat, Jan 31, 2009 at 12:12 PM, Andrei Alexandrescu
>> <SeeWebsiteForEmail at erdani.org> wrote:
>>>
>>> I've updated my code and documentation to include series (as in math) in
>>> the
>>> form of infinite ranges. Also series in closed form (given n can compute
>>> the
>>> nth value without iterating) are supported as random-access ranges.
>>>
>>> Also Stride is provided. The Matrix container (speaking of scientific
>>> computing with D!) will support various representational choices, most
>>> importantly the ones endorsed by high-performance libraries. For Matrix,
>>> Stride is an important component as I'm sure anyone who's ever written a
>>> matrix knows.
>>
>> This sounds interesting!  So you'll have diagonal, banded, upper/lower
>> triangular, symmetric, strided and dense matrices?  How about 2D
>> slices?
>
> Sparse matrices too. For dense block matrices I'm thinking of a type like
> this:

Whoa... excellent.  I didn't even dare ask about sparse, or block
sparse, or n-dimensional arrays.

> alias BlockMatrix!(float, 3, [0, 2, 1]) M3;
>
> This defines a block rectangular matrix with three dimensions laid out in
> the order specified: In this case, M3 is laid out as a collection of
> column-order submatrices of 2 dimensions. (I deliberately chose a slightly
> odd example.) The last parameter of BlockMatrix can be any permutation of 0,
> 1, and 2, thus dictating how the dimensions are laid out. For an obvious
> example, here's a Fortran-compatible block matrix with column-major layout:
>
> alias BlockMatrix!(double, 2, [1, 0]) FortranMatrix;

Nice.  Maybe you can have some predefined symbols for Fortran and C
layouts too, so users don't have to think about it.

> So it's pretty easy to choose any layout for any matrix with any number of
> dimensions. This matrix defines ranges that crawl it, called hyperplanes.
> There is a "natural" hyperplane that doesn't need any stride, and that's the
> cut along the dimension that varies slowest. That cut is the most efficient
> because it really has the layout of a (littler) BlockMatrix itself.
> Continuing on the M3 example, a "row" in that matrix has type:
>
> alias BlockMatrix!(float, 2, [1, 0]) M3Row;
>
> which is obtained by chopping the first element in the array [0, 2, 1] into
> [2, 1], and then decrementing what's left into [1, 2]. So if you have:
>
> M3 m;
> auto subM = m[3];
>
> you get that smaller BlockMatrix with all amenities its larger cousin
> offered. But if you want a cut along a different dimension, a strided range
> will be returned.
>
> Jagged, banded, diagonal, fixed-size, and sparse layouts come naturally and
> will be hopefully supported in a mixable-and-matchable form (e.g. I want a
> block matrix of 3x3 matrices). I want to focus on storage for now, put it in
> the standard library, to then allow great scientific programmers like you
> guys to use those storages as a common data format on top of which cool
> algorithms can be implemented.

Excellent plan.  This sounds like a great rallying point for numerics in D2.
But get those ranges done first!

--bb