Sorry for the follow-up, I just wanted to mention that I also uploaded
the proposal on Melange at
http://www.google-melange.com/gsoc/proposal/review/google/gsoc2012/cristicbz/41002#,
but I still recommend reading it on github.

Thank you,
Cristian.

---
Cristi Cobzarenco
BSc in Artificial Intelligence and Computer Science
University of Edinburgh
Profile: http://www.google.com/profiles/cristi.cobzarenco



On 27 March 2012 11:44, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
> Hello everyone,
>
> My name is Cristian Cobzarenco and last year, mentored by David Simcha and co-mentored by Fawzi Mohamed and Andrei Alexandrescu, I worked on a fork of SciD as part of Google's Summer of Code. While we've got a lot done last summer there's still plenty to do and my proposal this year is to get the library into a finished state where it can be reviewed for inclusion in Phobos as std.linalg. The easiest place to read my proposal is at the library's wiki: https://github.com/cristicbz/scid/wiki/GSoC-2012-Proposal, I will convert the formatting and add it to Melange soon.
>
> Feedback would be most welcome, I'm interested particularly in if you think the priorities I set are correct or if there is something you would (or wouldn't) rather see in the library. Let me know if anyone wants me to attach the content to an e-mail (it would be a pretty long e-mail), to allow in-line commenting.
>
> Yours,
> Cristian.
>
> ---
> Cristi Cobzarenco
> BSc in Artificial Intelligence and Computer Science
> University of Edinburgh
> Profile: http://www.google.com/profiles/cristi.cobzarenco

Yet another follow-up, sorry. I changed the proposal quite a bit
trying to do a better job at separating the very specific technical
content from the rest of the proposal which should make is easier to
read for people who didn't work on the library - this, I imagine, was
a hurdle for people who might have otherwise wanted to give feedback.
Also, in the meantime I found out there was a template for proposals,
having not seen it, I hadn't used. So the new structure reflects that
template more closely, too. Here's the link:
http://www.google-melange.com/gsoc/proposal/review/google/gsoc2012/cristicbz/41002#

---
Cristi Cobzarenco
BSc in Artificial Intelligence and Computer Science
University of Edinburgh
Profile: http://www.google.com/profiles/cristi.cobzarenco



On 27 March 2012 23:54, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
> Sorry for the follow-up, I just wanted to mention that I also uploaded
> the proposal on Melange at
> http://www.google-melange.com/gsoc/proposal/review/google/gsoc2012/cristicbz/41002#,
> but I still recommend reading it on github.
>
> Thank you,
> Cristian.
>
> ---
> Cristi Cobzarenco
> BSc in Artificial Intelligence and Computer Science
> University of Edinburgh
> Profile: http://www.google.com/profiles/cristi.cobzarenco
>
>
>
> On 27 March 2012 11:44, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
>> Hello everyone,
>>
>> My name is Cristian Cobzarenco and last year, mentored by David Simcha and co-mentored by Fawzi Mohamed and Andrei Alexandrescu, I worked on a fork of SciD as part of Google's Summer of Code. While we've got a lot done last summer there's still plenty to do and my proposal this year is to get the library into a finished state where it can be reviewed for inclusion in Phobos as std.linalg. The easiest place to read my proposal is at the library's wiki: https://github.com/cristicbz/scid/wiki/GSoC-2012-Proposal, I will convert the formatting and add it to Melange soon.
>>
>> Feedback would be most welcome, I'm interested particularly in if you think the priorities I set are correct or if there is something you would (or wouldn't) rather see in the library. Let me know if anyone wants me to attach the content to an e-mail (it would be a pretty long e-mail), to allow in-line commenting.
>>
>> Yours,
>> Cristian.
>>
>> ---
>> Cristi Cobzarenco
>> BSc in Artificial Intelligence and Computer Science
>> University of Edinburgh
>> Profile: http://www.google.com/profiles/cristi.cobzarenco

I've read **Proposed Changes and Additions**, and I would like to comment and ask a few questions if that's okay.  BTW, I've used Eigen a lot and I see some similarities here, but a direct rewrite may not be the best thing because D > C++.

2.  Change the matrix & vector types, adding fixed-sized matrix support in the process.

This is a step in the right direction I think, and by that I'm talking about the decision to remove the difference between a Vector and a Matrix.  Also, fixed-size matrices are also a must.  There is compile-time optimization that you won't be able to do for dynamic-size matrices.


3. Add value arrays (or numeric arrays, we can come up with a good name).

I really don't see the point for these.  We have the built-in arrays and the one in Phobos (which will get even better soon).


4. Add reductions, partial reductions, and broadcasting for matrices and arrays.

This one is similar to what we have in Eigen, but I don't understand why the operations are member functions (even in Eigen).  I much rather have something like this:

rowwise!sum(mat);

Also, that way the users can use their own custom functions with much ease.


6. Add support for interoperation with D built-in arrays (or pointers).

So I take that Matrix is not a sub-type?  why? If we have something like this:

struct Matrix(Real, size_t row, size_t col) {

  Real[row*col] data;
  alias data this;
}

then we wouldn't need any kind of interoperation with built-in arrays, would we?  I think this would save us a lot of headache.

That's just me and I could be wrong.

I've got to tell you though, I'm very excited about this project and I'll be watching it closely.

cheers.

On 3 April 2012 02:37, Caligo <iteronvexor@gmail.com> wrote:

> I've read **Proposed Changes and Additions**, and I would like to comment and ask a few questions if that's okay.  BTW, I've used Eigen a lot and I see some similarities here, but a direct rewrite may not be the best thing because D > C++.
>
> 2.  Change the matrix & vector types, adding fixed-sized matrix support in the process.
>
> This is a step in the right direction I think, and by that I'm talking about the decision to remove the difference between a Vector and a Matrix.  Also, fixed-size matrices are also a must.  There is compile-time optimization that you won't be able to do for dynamic-size matrices.
>
>
> 3. Add value arrays (or numeric arrays, we can come up with a good name).
>
> I really don't see the point for these.  We have the built-in arrays and the one in Phobos (which will get even better soon).
>

The point of these is to have light-weight element wise operation support.
It's true that in theory the built-in arrays do this. However, this library
is built on top BLAS/LAPACK, which means operations on large matrices will
be faster than on D arrays. Also, as far as I know, D doesn't support
allocating dynamic 2-D arrays (as in not arrays of arrays), not to mention
2-D slicing (keeping track of leading dimension).
Also I'm not sure how a case like this will be compiled, it may or may not
allocate a temporary:

a[] = b[] * c[] + d[] * 2.0;

The expression templates in SciD mean there will be no temporary allocation in this call.


>
> 4. Add reductions, partial reductions, and broadcasting for matrices and arrays.
>
> This one is similar to what we have in Eigen, but I don't understand why the operations are member functions (even in Eigen).  I much rather have something like this:
>
> rowwise!sum(mat);
>
> Also, that way the users can use their own custom functions with much ease.
>

There is a problem with this design. You want each matrix type (be it general, triangular, sparse or even an expression node)  do be able to define its own implementation of sum: calling the right BLAS function and making whatever specific optimisations they can. Since D doesn't have argument dependent look-up (ADL), users can't provide specialisations for their own types. The same arguments apply to rowwise() and columnwise() which will return proxies specific to the matrix type. You could do something like this, in principle:

auto sum( T )( T mat ) {
  return mat.sum();
}

And if we want that we can add it, but this will provide no addition in extensibility. By the way, you can use std.algorithm with matrices since they offer range functionality, but it will be much slower to use reduce!mySumFunction(mat) than mat.sum() which uses a BLAS backend.



>
>
> 6. Add support for interoperation with D built-in arrays (or pointers).
>
> So I take that Matrix is not a sub-type?  why? If we have something like this:
>
> struct Matrix(Real, size_t row, size_t col) {
>
>  Real[row*col] data;
>  alias data this;
> }
>
> then we wouldn't need any kind of interoperation with built-in arrays, would we?  I think this would save us a lot of headache.
>
> That's just me and I could be wrong.
>

Inter-operation referred more to having a matrix object wrapping a pointer to an already available piece of memory - maybe allocated through a region allocator, maybe resulting from some other library. This means we need to take care of different strides and different storage orders which cannot be handled by built-in arrays. Right now, matrices wrap ref-counted copy-on-write array types (ArrayData in the current code) - we decided last year that we don't want to use the garbage collector, because of its current issues. Also I would prefer not using the same Matrix type for pointer-wrappers and normal matrices because the former must have reference semantics while the latter have value semantics. I think it would be confusing if some matrices would copy their data and some would share memory.


>
> I've got to tell you though, I'm very excited about this project and I'll be watching it closely.
>
> cheers.
>
---
Cristi Cobzarenco
BSc in Artificial Intelligence and Computer Science
University of Edinburgh
Profile: http://www.google.com/profiles/cristi.cobzarenco

For the Point 4, I really like to have high order functions like
reduceRow and reduceCol. then the function argument is simply the
reduceRow!foo(0,mat), here the foo is not a function operating on the
whole column but simply a function of two elements (e.g.
reduceRow!("a+b")(0,mat)). Or even better we could have a general
reduce function with row and col as template parameter so that we can
do reduce!(foo,row)(0,mat). I dont know whether we can optimize such
reduce function for different matrix type, but such a function will be
extremely useful from a user perspective.

On Tue, Apr 3, 2012 at 7:20 PM, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
>
>
> On 3 April 2012 02:37, Caligo <iteronvexor@gmail.com> wrote:
>>
>> I've read **Proposed Changes and Additions**, and I would like to comment and ask a few questions if that's okay.  BTW, I've used Eigen a lot and I see some similarities here, but a direct rewrite may not be the best thing because D > C++.
>>
>> 2.  Change the matrix & vector types, adding fixed-sized matrix support in the process.
>>
>> This is a step in the right direction I think, and by that I'm talking about the decision to remove the difference between a Vector and a Matrix.  Also, fixed-size matrices are also a must.  There is compile-time optimization that you won't be able to do for dynamic-size matrices.
>>
>>
>> 3. Add value arrays (or numeric arrays, we can come up with a good name).
>>
>> I really don't see the point for these.  We have the built-in arrays and the one in Phobos (which will get even better soon).
>
>
> The point of these is to have light-weight element wise operation support.
> It's true that in theory the built-in arrays do this. However, this library
> is built on top BLAS/LAPACK, which means operations on large matrices will
> be faster than on D arrays. Also, as far as I know, D doesn't support
> allocating dynamic 2-D arrays (as in not arrays of arrays), not to mention
> 2-D slicing (keeping track of leading dimension).
> Also I'm not sure how a case like this will be compiled, it may or may not
> allocate a temporary:
>
> a[] = b[] * c[] + d[] * 2.0;
>
> The expression templates in SciD mean there will be no temporary allocation in this call.
>
>>
>>
>> 4. Add reductions, partial reductions, and broadcasting for matrices and arrays.
>>
>> This one is similar to what we have in Eigen, but I don't understand why the operations are member functions (even in Eigen).  I much rather have something like this:
>>
>> rowwise!sum(mat);
>>
>> Also, that way the users can use their own custom functions with much ease.
>
>
> There is a problem with this design. You want each matrix type (be it general, triangular, sparse or even an expression node)  do be able to define its own implementation of sum: calling the right BLAS function and making whatever specific optimisations they can. Since D doesn't have argument dependent look-up (ADL), users can't provide specialisations for their own types. The same arguments apply to rowwise() and columnwise() which will return proxies specific to the matrix type. You could do something like this, in principle:
>
> auto sum( T )( T mat ) {
>   return mat.sum();
> }
>
> And if we want that we can add it, but this will provide no addition in extensibility. By the way, you can use std.algorithm with matrices since they offer range functionality, but it will be much slower to use reduce!mySumFunction(mat) than mat.sum() which uses a BLAS backend.
>
>
>>
>>
>>
>> 6. Add support for interoperation with D built-in arrays (or pointers).
>>
>> So I take that Matrix is not a sub-type?  why? If we have something like this:
>>
>> struct Matrix(Real, size_t row, size_t col) {
>>
>>  Real[row*col] data;
>>  alias data this;
>> }
>>
>> then we wouldn't need any kind of interoperation with built-in arrays, would we?  I think this would save us a lot of headache.
>>
>> That's just me and I could be wrong.
>
>
> Inter-operation referred more to having a matrix object wrapping a pointer to an already available piece of memory - maybe allocated through a region allocator, maybe resulting from some other library. This means we need to take care of different strides and different storage orders which cannot be handled by built-in arrays. Right now, matrices wrap ref-counted copy-on-write array types (ArrayData in the current code) - we decided last year that we don't want to use the garbage collector, because of its current issues. Also I would prefer not using the same Matrix type for pointer-wrappers and normal matrices because the former must have reference semantics while the latter have value semantics. I think it would be confusing if some matrices would copy their data and some would share memory.
>
>>
>>
>> I've got to tell you though, I'm very excited about this project and I'll be watching it closely.
>>
>> cheers.
>
> ---
> Cristi Cobzarenco
> BSc in Artificial Intelligence and Computer Science
> University of Edinburgh
> Profile: http://www.google.com/profiles/cristi.cobzarenco
>

Btw, I really don't like the matrix api to be member functions. It is hard for user to extend the library in an unified way. It is also ugly when you want to chain the function calls.

On Wed, Apr 4, 2012 at 5:09 PM, Michael Chen <sth4nth@gmail.com> wrote:
> For the Point 4, I really like to have high order functions like
> reduceRow and reduceCol. then the function argument is simply the
> reduceRow!foo(0,mat), here the foo is not a function operating on the
> whole column but simply a function of two elements (e.g.
> reduceRow!("a+b")(0,mat)). Or even better we could have a general
> reduce function with row and col as template parameter so that we can
> do reduce!(foo,row)(0,mat). I dont know whether we can optimize such
> reduce function for different matrix type, but such a function will be
> extremely useful from a user perspective.
>
> On Tue, Apr 3, 2012 at 7:20 PM, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
>>
>>
>> On 3 April 2012 02:37, Caligo <iteronvexor@gmail.com> wrote:
>>>
>>> I've read **Proposed Changes and Additions**, and I would like to comment and ask a few questions if that's okay.  BTW, I've used Eigen a lot and I see some similarities here, but a direct rewrite may not be the best thing because D > C++.
>>>
>>> 2.  Change the matrix & vector types, adding fixed-sized matrix support in the process.
>>>
>>> This is a step in the right direction I think, and by that I'm talking about the decision to remove the difference between a Vector and a Matrix.  Also, fixed-size matrices are also a must.  There is compile-time optimization that you won't be able to do for dynamic-size matrices.
>>>
>>>
>>> 3. Add value arrays (or numeric arrays, we can come up with a good name).
>>>
>>> I really don't see the point for these.  We have the built-in arrays and the one in Phobos (which will get even better soon).
>>
>>
>> The point of these is to have light-weight element wise operation support.
>> It's true that in theory the built-in arrays do this. However, this library
>> is built on top BLAS/LAPACK, which means operations on large matrices will
>> be faster than on D arrays. Also, as far as I know, D doesn't support
>> allocating dynamic 2-D arrays (as in not arrays of arrays), not to mention
>> 2-D slicing (keeping track of leading dimension).
>> Also I'm not sure how a case like this will be compiled, it may or may not
>> allocate a temporary:
>>
>> a[] = b[] * c[] + d[] * 2.0;
>>
>> The expression templates in SciD mean there will be no temporary allocation in this call.
>>
>>>
>>>
>>> 4. Add reductions, partial reductions, and broadcasting for matrices and arrays.
>>>
>>> This one is similar to what we have in Eigen, but I don't understand why the operations are member functions (even in Eigen).  I much rather have something like this:
>>>
>>> rowwise!sum(mat);
>>>
>>> Also, that way the users can use their own custom functions with much ease.
>>
>>
>> There is a problem with this design. You want each matrix type (be it general, triangular, sparse or even an expression node)  do be able to define its own implementation of sum: calling the right BLAS function and making whatever specific optimisations they can. Since D doesn't have argument dependent look-up (ADL), users can't provide specialisations for their own types. The same arguments apply to rowwise() and columnwise() which will return proxies specific to the matrix type. You could do something like this, in principle:
>>
>> auto sum( T )( T mat ) {
>>   return mat.sum();
>> }
>>
>> And if we want that we can add it, but this will provide no addition in extensibility. By the way, you can use std.algorithm with matrices since they offer range functionality, but it will be much slower to use reduce!mySumFunction(mat) than mat.sum() which uses a BLAS backend.
>>
>>
>>>
>>>
>>>
>>> 6. Add support for interoperation with D built-in arrays (or pointers).
>>>
>>> So I take that Matrix is not a sub-type?  why? If we have something like this:
>>>
>>> struct Matrix(Real, size_t row, size_t col) {
>>>
>>>  Real[row*col] data;
>>>  alias data this;
>>> }
>>>
>>> then we wouldn't need any kind of interoperation with built-in arrays, would we?  I think this would save us a lot of headache.
>>>
>>> That's just me and I could be wrong.
>>
>>
>> Inter-operation referred more to having a matrix object wrapping a pointer to an already available piece of memory - maybe allocated through a region allocator, maybe resulting from some other library. This means we need to take care of different strides and different storage orders which cannot be handled by built-in arrays. Right now, matrices wrap ref-counted copy-on-write array types (ArrayData in the current code) - we decided last year that we don't want to use the garbage collector, because of its current issues. Also I would prefer not using the same Matrix type for pointer-wrappers and normal matrices because the former must have reference semantics while the latter have value semantics. I think it would be confusing if some matrices would copy their data and some would share memory.
>>
>>>
>>>
>>> I've got to tell you though, I'm very excited about this project and I'll be watching it closely.
>>>
>>> cheers.
>>
>> ---
>> Cristi Cobzarenco
>> BSc in Artificial Intelligence and Computer Science
>> University of Edinburgh
>> Profile: http://www.google.com/profiles/cristi.cobzarenco
>>

another btw, there is also another great c++ linear algebra library
besides Eigen: Amardillo which has very simple matlab like interface
 and performs on a par with Eigen.
On Wed, Apr 4, 2012 at 5:14 PM, Michael Chen <sth4nth@gmail.com> wrote:
> Btw, I really don't like the matrix api to be member functions. It is hard for user to extend the library in an unified way. It is also ugly when you want to chain the function calls.
>
> On Wed, Apr 4, 2012 at 5:09 PM, Michael Chen <sth4nth@gmail.com> wrote:
>> For the Point 4, I really like to have high order functions like
>> reduceRow and reduceCol. then the function argument is simply the
>> reduceRow!foo(0,mat), here the foo is not a function operating on the
>> whole column but simply a function of two elements (e.g.
>> reduceRow!("a+b")(0,mat)). Or even better we could have a general
>> reduce function with row and col as template parameter so that we can
>> do reduce!(foo,row)(0,mat). I dont know whether we can optimize such
>> reduce function for different matrix type, but such a function will be
>> extremely useful from a user perspective.
>>
>> On Tue, Apr 3, 2012 at 7:20 PM, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
>>>
>>>
>>> On 3 April 2012 02:37, Caligo <iteronvexor@gmail.com> wrote:
>>>>
>>>> I've read **Proposed Changes and Additions**, and I would like to comment and ask a few questions if that's okay.  BTW, I've used Eigen a lot and I see some similarities here, but a direct rewrite may not be the best thing because D > C++.
>>>>
>>>> 2.  Change the matrix & vector types, adding fixed-sized matrix support in the process.
>>>>
>>>> This is a step in the right direction I think, and by that I'm talking about the decision to remove the difference between a Vector and a Matrix.  Also, fixed-size matrices are also a must.  There is compile-time optimization that you won't be able to do for dynamic-size matrices.
>>>>
>>>>
>>>> 3. Add value arrays (or numeric arrays, we can come up with a good name).
>>>>
>>>> I really don't see the point for these.  We have the built-in arrays and the one in Phobos (which will get even better soon).
>>>
>>>
>>> The point of these is to have light-weight element wise operation support.
>>> It's true that in theory the built-in arrays do this. However, this library
>>> is built on top BLAS/LAPACK, which means operations on large matrices will
>>> be faster than on D arrays. Also, as far as I know, D doesn't support
>>> allocating dynamic 2-D arrays (as in not arrays of arrays), not to mention
>>> 2-D slicing (keeping track of leading dimension).
>>> Also I'm not sure how a case like this will be compiled, it may or may not
>>> allocate a temporary:
>>>
>>> a[] = b[] * c[] + d[] * 2.0;
>>>
>>> The expression templates in SciD mean there will be no temporary allocation in this call.
>>>
>>>>
>>>>
>>>> 4. Add reductions, partial reductions, and broadcasting for matrices and arrays.
>>>>
>>>> This one is similar to what we have in Eigen, but I don't understand why the operations are member functions (even in Eigen).  I much rather have something like this:
>>>>
>>>> rowwise!sum(mat);
>>>>
>>>> Also, that way the users can use their own custom functions with much ease.
>>>
>>>
>>> There is a problem with this design. You want each matrix type (be it general, triangular, sparse or even an expression node)  do be able to define its own implementation of sum: calling the right BLAS function and making whatever specific optimisations they can. Since D doesn't have argument dependent look-up (ADL), users can't provide specialisations for their own types. The same arguments apply to rowwise() and columnwise() which will return proxies specific to the matrix type. You could do something like this, in principle:
>>>
>>> auto sum( T )( T mat ) {
>>>   return mat.sum();
>>> }
>>>
>>> And if we want that we can add it, but this will provide no addition in extensibility. By the way, you can use std.algorithm with matrices since they offer range functionality, but it will be much slower to use reduce!mySumFunction(mat) than mat.sum() which uses a BLAS backend.
>>>
>>>
>>>>
>>>>
>>>>
>>>> 6. Add support for interoperation with D built-in arrays (or pointers).
>>>>
>>>> So I take that Matrix is not a sub-type?  why? If we have something like this:
>>>>
>>>> struct Matrix(Real, size_t row, size_t col) {
>>>>
>>>>  Real[row*col] data;
>>>>  alias data this;
>>>> }
>>>>
>>>> then we wouldn't need any kind of interoperation with built-in arrays, would we?  I think this would save us a lot of headache.
>>>>
>>>> That's just me and I could be wrong.
>>>
>>>
>>> Inter-operation referred more to having a matrix object wrapping a pointer to an already available piece of memory - maybe allocated through a region allocator, maybe resulting from some other library. This means we need to take care of different strides and different storage orders which cannot be handled by built-in arrays. Right now, matrices wrap ref-counted copy-on-write array types (ArrayData in the current code) - we decided last year that we don't want to use the garbage collector, because of its current issues. Also I would prefer not using the same Matrix type for pointer-wrappers and normal matrices because the former must have reference semantics while the latter have value semantics. I think it would be confusing if some matrices would copy their data and some would share memory.
>>>
>>>>
>>>>
>>>> I've got to tell you though, I'm very excited about this project and I'll be watching it closely.
>>>>
>>>> cheers.
>>>
>>> ---
>>> Cristi Cobzarenco
>>> BSc in Artificial Intelligence and Computer Science
>>> University of Edinburgh
>>> Profile: http://www.google.com/profiles/cristi.cobzarenco
>>>

Thanks for the feedback!

On 4 April 2012 10:21, Michael Chen <sth4nth@gmail.com> wrote:

> another btw, there is also another great c++ linear algebra library
> besides Eigen: Amardillo which has very simple matlab like interface
>  and performs on a par with Eigen.
>

I'll look into it, thanks.


> On Wed, Apr 4, 2012 at 5:14 PM, Michael Chen <sth4nth@gmail.com> wrote:
> > Btw, I really don't like the matrix api to be member functions. It is hard for user to extend the library in an unified way. It is also ugly when you want to chain the function calls.
>
>
> > On Wed, Apr 4, 2012 at 5:09 PM, Michael Chen <sth4nth@gmail.com> wrote:
> >> For the Point 4, I really like to have high order functions like
> >> reduceRow and reduceCol. then the function argument is simply the
> >> reduceRow!foo(0,mat), here the foo is not a function operating on the
> >> whole column but simply a function of two elements (e.g.
> >> reduceRow!("a+b")(0,mat)). Or even better we could have a general
> >> reduce function with row and col as template parameter so that we can
> >> do reduce!(foo,row)(0,mat). I dont know whether we can optimize such
> >> reduce function for different matrix type, but such a function will be
> >> extremely useful from a user perspective
>

Well, as I said before, there's nothing stopping us from providing the free functions that call the member functionsl. However there is something stopping us from not providing a member function alternative. D's lack of ADL means we would not allow easy extensibility of possible operations. Say Alice invents a new kind of matrix type, the DiagonalMatrix type which stores its elements in a 1D array (this isn't exactly how it would work with the design we have, she would in fact have to define a storage type, but bear with me). If she wants to implement sum on her matrix, she can't simply add a specialisation to the sum( T ) function because of the lack of ADL. If instead we implemented sum(T) as:

auto sum( T )( T matrix ) {
     static if( is( typeof( T.init.sum() ) ) )
         return matrix.sum;
     else
         return reduce!"a+b"( 0, matrix.elements() );
}

then Alice could simply define a DiagonalMatrix.sum() method that wouldn't need to go through all the zero elements, for example. The idea of rowReduce and columnReduce still works - we can provide this for when a user wants to use her own operation for reducing. We could also have a way of automatically optimising rowReduce!"a+b"( mat ) by calling mat.rowwise().sum() if the operation is available - but that wouldn't be exactly high priority.


> .
> >>
> >> On Tue, Apr 3, 2012 at 7:20 PM, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
> >>>
> >>>
> >>> On 3 April 2012 02:37, Caligo <iteronvexor@gmail.com> wrote:
> >>>>
> >>>> I've read **Proposed Changes and Additions**, and I would like to comment and ask a few questions if that's okay.  BTW, I've used Eigen a lot and I see some similarities here, but a direct rewrite may not be the best thing because D > C++.
> >>>>
> >>>> 2.  Change the matrix & vector types, adding fixed-sized matrix support in the process.
> >>>>
> >>>> This is a step in the right direction I think, and by that I'm talking about the decision to remove the difference between a Vector and a Matrix.  Also, fixed-size matrices are also a must.  There is compile-time optimization that you won't be able to do for dynamic-size matrices.
> >>>>
> >>>>
> >>>> 3. Add value arrays (or numeric arrays, we can come up with a good
> name).
> >>>>
> >>>> I really don't see the point for these.  We have the built-in arrays and the one in Phobos (which will get even better soon).
> >>>
> >>>
> >>> The point of these is to have light-weight element wise operation
> support.
> >>> It's true that in theory the built-in arrays do this. However, this
> library
> >>> is built on top BLAS/LAPACK, which means operations on large matrices
> will
> >>> be faster than on D arrays. Also, as far as I know, D doesn't support allocating dynamic 2-D arrays (as in not arrays of arrays), not to
> mention
> >>> 2-D slicing (keeping track of leading dimension).
> >>> Also I'm not sure how a case like this will be compiled, it may or may
> not
> >>> allocate a temporary:
> >>>
> >>> a[] = b[] * c[] + d[] * 2.0;
> >>>
> >>> The expression templates in SciD mean there will be no temporary
> allocation
> >>> in this call.
> >>>
> >>>>
> >>>>
> >>>> 4. Add reductions, partial reductions, and broadcasting for matrices
> and
> >>>> arrays.
> >>>>
> >>>> This one is similar to what we have in Eigen, but I don't understand why the operations are member functions (even in Eigen).  I much rather have something like this:
> >>>>
> >>>> rowwise!sum(mat);
> >>>>
> >>>> Also, that way the users can use their own custom functions with much ease.
> >>>
> >>>
> >>> There is a problem with this design. You want each matrix type (be it general, triangular, sparse or even an expression node)  do be able to define its own implementation of sum: calling the right BLAS function
> and
> >>> making whatever specific optimisations they can. Since D doesn't have argument dependent look-up (ADL), users can't provide specialisations
> for
> >>> their own types. The same arguments apply to rowwise() and columnwise()
> >>> which will return proxies specific to the matrix type. You could do
> >>> something like this, in principle:
> >>>
> >>> auto sum( T )( T mat ) {
> >>>   return mat.sum();
> >>> }
> >>>
> >>> And if we want that we can add it, but this will provide no addition in extensibility. By the way, you can use std.algorithm with matrices
> since
> >>> they offer range functionality, but it will be much slower to use
> >>> reduce!mySumFunction(mat) than mat.sum() which uses a BLAS backend.
> >>>
> >>>
> >>>>
> >>>>
> >>>>
> >>>> 6. Add support for interoperation with D built-in arrays (or
> pointers).
> >>>>
> >>>> So I take that Matrix is not a sub-type?  why? If we have something
> like
> >>>> this:
> >>>>
> >>>> struct Matrix(Real, size_t row, size_t col) {
> >>>>
> >>>>  Real[row*col] data;
> >>>>  alias data this;
> >>>> }
> >>>>
> >>>> then we wouldn't need any kind of interoperation with built-in arrays, would we?  I think this would save us a lot of headache.
> >>>>
> >>>> That's just me and I could be wrong.
> >>>
> >>>
> >>> Inter-operation referred more to having a matrix object wrapping a
> pointer
> >>> to an already available piece of memory - maybe allocated through a
> region
> >>> allocator, maybe resulting from some other library. This means we need
> to
> >>> take care of different strides and different storage orders which
> cannot be
> >>> handled by built-in arrays. Right now, matrices wrap ref-counted copy-on-write array types (ArrayData in the current code) - we decided
> last
> >>> year that we don't want to use the garbage collector, because of its
> current
> >>> issues. Also I would prefer not using the same Matrix type for pointer-wrappers and normal matrices because the former must have
> reference
> >>> semantics while the latter have value semantics. I think it would be confusing if some matrices would copy their data and some would share memory.
> >>>
> >>>>
> >>>>
> >>>> I've got to tell you though, I'm very excited about this project and I'll be watching it closely.
> >>>>
> >>>> cheers.
> >>>
> >>> ---
> >>> Cristi Cobzarenco
> >>> BSc in Artificial Intelligence and Computer Science
> >>> University of Edinburgh
> >>> Profile: http://www.google.com/profiles/cristi.cobzarenco
> >>>
>

Thanks for the explanation, now I get it. In case you are interested,
there is excellent article about monad style c++ template meta
programming by Bartosz Milewski which might be helpful for compile
time optimization for evaluation order.
Really looking forward to official release of the SciD. D is really
suitable for scientific computation. It will be great to have an
efficient and easy-to-use linear algebra library.


On Thu, Apr 5, 2012 at 7:42 AM, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
> Thanks for the feedback!
>
> On 4 April 2012 10:21, Michael Chen <sth4nth@gmail.com> wrote:
>>
>> another btw, there is also another great c++ linear algebra library besides Eigen: Amardillo which has very simple matlab like interface  and performs on a par with Eigen.
>
>
> I'll look into it, thanks.
>
>>
>> On Wed, Apr 4, 2012 at 5:14 PM, Michael Chen <sth4nth@gmail.com> wrote:
>> > Btw, I really don't like the matrix api to be member functions. It is hard for user to extend the library in an unified way. It is also ugly when you want to chain the function calls.
>>
>> >
>> > On Wed, Apr 4, 2012 at 5:09 PM, Michael Chen <sth4nth@gmail.com> wrote:
>> >> For the Point 4, I really like to have high order functions like
>> >> reduceRow and reduceCol. then the function argument is simply the
>> >> reduceRow!foo(0,mat), here the foo is not a function operating on the
>> >> whole column but simply a function of two elements (e.g.
>> >> reduceRow!("a+b")(0,mat)). Or even better we could have a general
>> >> reduce function with row and col as template parameter so that we can
>> >> do reduce!(foo,row)(0,mat). I dont know whether we can optimize such
>> >> reduce function for different matrix type, but such a function will be
>> >> extremely useful from a user perspective
>
>
> Well, as I said before, there's nothing stopping us from providing the free functions that call the member functionsl. However there is something stopping us from not providing a member function alternative. D's lack of ADL means we would not allow easy extensibility of possible operations. Say Alice invents a new kind of matrix type, the DiagonalMatrix type which stores its elements in a 1D array (this isn't exactly how it would work with the design we have, she would in fact have to define a storage type, but bear with me). If she wants to implement sum on her matrix, she can't simply add a specialisation to the sum( T ) function because of the lack of ADL. If instead we implemented sum(T) as:
>
> auto sum( T )( T matrix ) {
>      static if( is( typeof( T.init.sum() ) ) )
>          return matrix.sum;
>      else
>          return reduce!"a+b"( 0, matrix.elements() );
> }
>
> then Alice could simply define a DiagonalMatrix.sum() method that wouldn't need to go through all the zero elements, for example. The idea of rowReduce and columnReduce still works - we can provide this for when a user wants to use her own operation for reducing. We could also have a way of automatically optimising rowReduce!"a+b"( mat ) by calling mat.rowwise().sum() if the operation is available - but that wouldn't be exactly high priority.
>
>>
>> .
>> >>
>> >> On Tue, Apr 3, 2012 at 7:20 PM, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
>> >>>
>> >>>
>> >>> On 3 April 2012 02:37, Caligo <iteronvexor@gmail.com> wrote:
>> >>>>
>> >>>> I've read **Proposed Changes and Additions**, and I would like to comment and ask a few questions if that's okay.  BTW, I've used Eigen a lot and I see some similarities here, but a direct rewrite may not be the best thing because D > C++.
>> >>>>
>> >>>> 2.  Change the matrix & vector types, adding fixed-sized matrix support in the process.
>> >>>>
>> >>>> This is a step in the right direction I think, and by that I'm
>> >>>> talking
>> >>>> about the decision to remove the difference between a Vector and a
>> >>>> Matrix.  Also, fixed-size matrices are also a must.  There is
>> >>>> compile-time optimization that you won't be able to do for
>> >>>> dynamic-size matrices.
>> >>>>
>> >>>>
>> >>>> 3. Add value arrays (or numeric arrays, we can come up with a good
>> >>>> name).
>> >>>>
>> >>>> I really don't see the point for these.  We have the built-in arrays and the one in Phobos (which will get even better soon).
>> >>>
>> >>>
>> >>> The point of these is to have light-weight element wise operation
>> >>> support.
>> >>> It's true that in theory the built-in arrays do this. However, this
>> >>> library
>> >>> is built on top BLAS/LAPACK, which means operations on large matrices
>> >>> will
>> >>> be faster than on D arrays. Also, as far as I know, D doesn't support
>> >>> allocating dynamic 2-D arrays (as in not arrays of arrays), not to
>> >>> mention
>> >>> 2-D slicing (keeping track of leading dimension).
>> >>> Also I'm not sure how a case like this will be compiled, it may or may
>> >>> not
>> >>> allocate a temporary:
>> >>>
>> >>> a[] = b[] * c[] + d[] * 2.0;
>> >>>
>> >>> The expression templates in SciD mean there will be no temporary
>> >>> allocation
>> >>> in this call.
>> >>>
>> >>>>
>> >>>>
>> >>>> 4. Add reductions, partial reductions, and broadcasting for matrices
>> >>>> and
>> >>>> arrays.
>> >>>>
>> >>>> This one is similar to what we have in Eigen, but I don't understand why the operations are member functions (even in Eigen).  I much rather have something like this:
>> >>>>
>> >>>> rowwise!sum(mat);
>> >>>>
>> >>>> Also, that way the users can use their own custom functions with much ease.
>> >>>
>> >>>
>> >>> There is a problem with this design. You want each matrix type (be it
>> >>> general, triangular, sparse or even an expression node)  do be able to
>> >>> define its own implementation of sum: calling the right BLAS function
>> >>> and
>> >>> making whatever specific optimisations they can. Since D doesn't have
>> >>> argument dependent look-up (ADL), users can't provide specialisations
>> >>> for
>> >>> their own types. The same arguments apply to rowwise() and
>> >>> columnwise()
>> >>> which will return proxies specific to the matrix type. You could do
>> >>> something like this, in principle:
>> >>>
>> >>> auto sum( T )( T mat ) {
>> >>>   return mat.sum();
>> >>> }
>> >>>
>> >>> And if we want that we can add it, but this will provide no addition
>> >>> in
>> >>> extensibility. By the way, you can use std.algorithm with matrices
>> >>> since
>> >>> they offer range functionality, but it will be much slower to use
>> >>> reduce!mySumFunction(mat) than mat.sum() which uses a BLAS backend.
>> >>>
>> >>>
>> >>>>
>> >>>>
>> >>>>
>> >>>> 6. Add support for interoperation with D built-in arrays (or
>> >>>> pointers).
>> >>>>
>> >>>> So I take that Matrix is not a sub-type?  why? If we have something
>> >>>> like
>> >>>> this:
>> >>>>
>> >>>> struct Matrix(Real, size_t row, size_t col) {
>> >>>>
>> >>>>  Real[row*col] data;
>> >>>>  alias data this;
>> >>>> }
>> >>>>
>> >>>> then we wouldn't need any kind of interoperation with built-in
>> >>>> arrays,
>> >>>> would we?  I think this would save us a lot of headache.
>> >>>>
>> >>>> That's just me and I could be wrong.
>> >>>
>> >>>
>> >>> Inter-operation referred more to having a matrix object wrapping a
>> >>> pointer
>> >>> to an already available piece of memory - maybe allocated through a
>> >>> region
>> >>> allocator, maybe resulting from some other library. This means we need
>> >>> to
>> >>> take care of different strides and different storage orders which
>> >>> cannot be
>> >>> handled by built-in arrays. Right now, matrices wrap ref-counted
>> >>> copy-on-write array types (ArrayData in the current code) - we decided
>> >>> last
>> >>> year that we don't want to use the garbage collector, because of its
>> >>> current
>> >>> issues. Also I would prefer not using the same Matrix type for
>> >>> pointer-wrappers and normal matrices because the former must have
>> >>> reference
>> >>> semantics while the latter have value semantics. I think it would be
>> >>> confusing if some matrices would copy their data and some would share
>> >>> memory.
>> >>>
>> >>>>
>> >>>>
>> >>>> I've got to tell you though, I'm very excited about this project and I'll be watching it closely.
>> >>>>
>> >>>> cheers.
>> >>>
>> >>> ---
>> >>> Cristi Cobzarenco
>> >>> BSc in Artificial Intelligence and Computer Science
>> >>> University of Edinburgh
>> >>> Profile: http://www.google.com/profiles/cristi.cobzarenco
>> >>>
>
>

On Tue, Apr 3, 2012 at 6:20 AM, Cristi Cobzarenco <cristi.cobzarenco@gmail.com> wrote:
>
> The point of these is to have light-weight element wise operation support. It's true that in theory the built-in arrays do this. However, this library is built on top BLAS/LAPACK, which means operations on large matrices will be faster than on D arrays.
>

I can't agree with building it on top of LAPACK or any other BLAS implementation, but perhaps I shouldn't complain because I'm not the one who's going to be implementing it.  I like the approach Eigen has taken where it offers its own BLAS implementation and, iirc, other BLAS libraries can be used as optional back-ends.

> Also, as far as I know, D doesn't support
> allocating dynamic 2-D arrays (as in not arrays of arrays), not to mention
> 2-D slicing (keeping track of leading dimension).
>

I fail to see why there is any need for 2D arrays.  We need to make sure multidimensional arrays (matrices) have data in very good arrangements.  This is called tiling and it requires 1D arrays:  2D arrays are stored as 1D arrays together with an indexing mechanism.

> Also I'm not sure how a case like this will be compiled, it may or may not allocate a temporary:
>
> a[] = b[] * c[] + d[] * 2.0;
>
> The expression templates in SciD mean there will be no temporary allocation in this call.
>

Why are expression templates used?  Are you pretty much rewriting Eigen in D ?  I don't see why expression templates would be needed with D.