I stumbled upon https://software.intel.com/en-us/node/471374 which gives good detail on Intel's Math Kernel Library's data formats for sparse matrices.

No doubt other popular linear algebra libraries have similar documentation. I was thinking we could start with adding these layouts to std, along with a few simple primitives (construction, element/slice access, stride etc). Then, people may just use those as they are or link with the linalg libraries for specific computations.


Thoughts?

Andrei

On 15/08/2015 2:57 a.m., Andrei Alexandrescu wrote:
> I stumbled upon https://software.intel.com/en-us/node/471374 which gives
> good detail on Intel's Math Kernel Library's data formats for sparse
> matrices.
>
> No doubt other popular linear algebra libraries have similar
> documentation. I was thinking we could start with adding these layouts
> to std, along with a few simple primitives (construction, element/slice
> access, stride etc). Then, people may just use those as they are or link
> with the linalg libraries for specific computations.
>
>
> Thoughts?
>
> Andrei

What would we gain other gl3n in Phobos?

On Friday, 14 August 2015 at 14:57:19 UTC, Andrei Alexandrescu wrote:
> I stumbled upon https://software.intel.com/en-us/node/471374 which gives good detail on Intel's Math Kernel Library's data formats for sparse matrices.
>
> No doubt other popular linear algebra libraries have similar documentation. I was thinking we could start with adding these layouts to std, along with a few simple primitives (construction, element/slice access, stride etc). Then, people may just use those as they are or link with the linalg libraries for specific computations.
>
>
> Thoughts?
>
> Andrei

Are sparse matrices a common scenario?

If anything I think small vectors, non-sparse matrixes and rectangles/AABB could be part of Phobos before that.

On 8/14/15 11:11 AM, ponce wrote:
> On Friday, 14 August 2015 at 14:57:19 UTC, Andrei Alexandrescu wrote:
>> I stumbled upon https://software.intel.com/en-us/node/471374 which
>> gives good detail on Intel's Math Kernel Library's data formats for
>> sparse matrices.
>>
>> No doubt other popular linear algebra libraries have similar
>> documentation. I was thinking we could start with adding these layouts
>> to std, along with a few simple primitives (construction,
>> element/slice access, stride etc). Then, people may just use those as
>> they are or link with the linalg libraries for specific computations.
>>
>>
>> Thoughts?
>>
>> Andrei
>
> Are sparse matrices a common scenario?

In machine learning, yes. But order is not as important as the basic strategy: layouts + simple/naive primitives in Phobos, elaborate primitives in specialized external libraries. -- Andrei

On Friday, 14 August 2015 at 15:11:39 UTC, ponce wrote:
> On Friday, 14 August 2015 at 14:57:19 UTC, Andrei Alexandrescu wrote:
>> I stumbled upon https://software.intel.com/en-us/node/471374 which gives good detail on Intel's Math Kernel Library's data formats for sparse matrices.
>>
>> No doubt other popular linear algebra libraries have similar documentation. I was thinking we could start with adding these layouts to std, along with a few simple primitives (construction, element/slice access, stride etc). Then, people may just use those as they are or link with the linalg libraries for specific computations.
>>
>>
>> Thoughts?
>>
>> Andrei
>
> Are sparse matrices a common scenario?
>
> If anything I think small vectors, non-sparse matrixes and rectangles/AABB could be part of Phobos before that.

+1

On Friday, 14 August 2015 at 15:11:39 UTC, ponce wrote:
> Are sparse matrices a common scenario?

Yes. They tend to pop up in virtually all "serious" numerical problems in science and engineering, particularly when partial differential equations are involved.

> If anything I think small vectors, non-sparse matrixes and rectangles/AABB could be part of Phobos before that.

If you just had a go at it from the CG/gamedev perspective, you'd probably end up with an API that is entirely unsuitable for larger problems. This might not be a big issue (just have two separate APIs), but we'd need to make it a conscious decision.

 — David

On Friday, 14 August 2015 at 18:51:51 UTC, David Nadlinger wrote:
> On Friday, 14 August 2015 at 15:11:39 UTC, ponce wrote:
>> Are sparse matrices a common scenario?
>
> Yes. They tend to pop up in virtually all "serious" numerical problems in science and engineering, particularly when partial differential equations are involved.
>
>> If anything I think small vectors, non-sparse matrixes and rectangles/AABB could be part of Phobos before that.
>
> If you just had a go at it from the CG/gamedev perspective, you'd probably end up with an API that is entirely unsuitable for larger problems. This might not be a big issue (just have two separate APIs), but we'd need to make it a conscious decision.
>
>  — David

Personally I would prefer small fixed sized vectors & matrices with game dev focused api be kept separate from a big/sparse matrix api.

On Friday, 14 August 2015 at 14:57:19 UTC, Andrei Alexandrescu wrote:
> I stumbled upon https://software.intel.com/en-us/node/471374 which gives good detail on Intel's Math Kernel Library's data formats for sparse matrices.
>
> No doubt other popular linear algebra libraries have similar documentation. I was thinking we could start with adding these layouts to std, along with a few simple primitives (construction, element/slice access, stride etc). Then, people may just use those as they are or link with the linalg libraries for specific computations.
>
>
> Thoughts?
>
> Andrei

One concern I have is the choice of MKL, which due to cost and license reasons, many developers will not have on all (or even any) of their machines.

I don't work with sparse matrices often so I do not know which libraries are most popular, but at a minimum I think it is necessary to say you can use it with a popular open source library. Given the importance of MKL, it would be a bad idea to not offer a compatible format, but it would be equally bad to focus on only MKL.

On Friday, 14 August 2015 at 14:57:19 UTC, Andrei Alexandrescu wrote:
> I stumbled upon https://software.intel.com/en-us/node/471374 which gives good detail on Intel's Math Kernel Library's data formats for sparse matrices.
>
> No doubt other popular linear algebra libraries have similar documentation. I was thinking we could start with adding these layouts to std, along with a few simple primitives (construction, element/slice access, stride etc). Then, people may just use those as they are or link with the linalg libraries for specific computations.
>
>
> Thoughts?
>
> Andrei

I am very excited that there is a greater focus on this.

I have been following this for a little while:
https://github.com/D-Programming-Language/phobos/pull/3397
It already does or intends to cover many things that I think are needed: easily looping through every element of a nd-slice, looping by dimension (by row or by column), and deep map. I wouldn't be surprised if there are language improvements that can facilitate this work.

In terms of formatting, I was playing around with a hybrid array struct (below) that is basically a static array, but also sort of range-ified. I was using the concept of domains from Chapel, albeit at an elementary level that could be significantly improved.

For that matter, may I suggest that you look over the Chapel specification:
http://chapel.cray.com/spec/spec-0.91.pdf
most importantly the domains and arrays sections. From writing that hybrid_array code, I think D lends itself well to some of these ideas.

Nevertheless, I agree that interoperability with C-based linear algebra libraries is important. If anything I said would make it more difficult, then ignore it.

import std.traits;
import std.range;
import std.array : array;
import std.stdio : writeln;
import std.algorithm : map;

struct hybrid_array(T : U[N], U, size_t N)
	if ( isStaticArray!T )
{
	T static_array;
	private auto _length = static_array.length;
	auto domain = iota(0, static_array.length);
	
	this(T x)
	{
		static_array = x;
	}
	
	@property bool empty()
    {
        return domain.empty;
    }

    @property size_t length()
    {
        return _length;
    }

    @property U front()
    {
		assert(!empty);
        return static_array[domain.front];
    }

    void popFront()
    {
		assert(!empty);
        domain.popFront;
		_length--;
    }
	
	@property hybrid_array save()
    {
        return this;
    }
	
	@property U back()
    {
		assert(!empty);
        return static_array[domain.back];
    }

    void popBack()
    {
		assert(!empty);
        domain.popBack;
		_length++;
    }
	
	U opIndex(size_t val)
    {
		assert(!empty);
        return static_array[domain[val]];
    }
}

void main()
{
	int[5] x = [0, 10, 2, 6, 1];
	//auto squares = map!(a => a * a)(x); //doesn't work
	
	auto range = iota(0, x.length).array;
	
	alias h_array = hybrid_array!(int[5]);
	auto y = h_array(x);
	
	writeln( isRandomAccessRange!(typeof(y)) );
	
	auto squares = map!(a => a * a)(y); //success, does work
	writeln(squares);
}

On Friday, 14 August 2015 at 18:51:51 UTC, David Nadlinger wrote:
>
>> If anything I think small vectors, non-sparse matrixes and rectangles/AABB could be part of Phobos before that.
>
> If you just had a go at it from the CG/gamedev perspective, you'd probably end up with an API that is entirely unsuitable for larger problems. This might not be a big issue (just have two separate APIs), but we'd need to make it a conscious decision.

While I'm interested in working with big (500,000+) dense (and sometimes smaller sparse) rectangular matrices, my sense of the D forum is that there are more CG/gamedev people using D.

Nevertheless, I think it is important to get the general structure right in a way that can be specialized to work with small matrices, rectangular matrices, strided matrices, or sparse matrices. I don't see an issue with prioritizing what users want in terms of fleshing out the APIs, so long as the different ways people use matrix libraries is in the back of the designers minds when they design it.

I was looking over some of the OpenGL matrix libraries recently. They have types like mat4 or vec3. I couldn't help but think that template constraints (at least for static arrays) would probably make it much easier. For instance, if you have a function that takes a vec, you can specialize it to do different things depending on if the array is small or large (small, unroll loops, large use parallelism, maybe with the decision on when stop unrolling loops or starting to use parallelism determined by AEOS, as in ATLAS).