Andrei Alexandrescu wrote:
Lars Kyllingstad wrote:
I am writing a D library based some of the stuff in SLATEC, and I've
come to a point where I need to decide on a way to manipulate vectors
and matrices. To that end, I have some ideas and questions I would
like comments on from the community.
Ideally, I want to restrict the user as little as possible, so I'm
writing heavily templated code in which one can use both
library-defined vector/matrix types and built-in arrays (both static
and dynamic). My reasons for this are:
a) Different problems may benefit from different types. Sparse
matrices, dense matrices, triangular matrices, etc. can all be
represented differently based on efficiency and/or memory requirements.
I use all of the above. It would be great to have them all within an
integrated framework.
I don't think I'm the man to provide you with that, at least not yet. I
have next to no experience with high-performance linear algebra. This is
a major part of the reason why I want to let people choose for
themselves what matrix types/libraries they want to use in conjunction
with my stuff.
Therefore I am, for now, focusing on other things. I am about halfway
done writing D versions of the QUADPACK routines, and I have started
work on MINPACK. It is for the latter that some basic linear algebra
functionality is needed.
You could take a look at Bill Baxter's MultiArray library, which
contains wrappers for several linear algebra libraries, as well as
storage formats for different types of matrices. Also, I think, it works
with Don's BLADE library. (Kudos to Don for that awesome name, by the
way -- it's a lot cooler than BLAS.) Both are written in D1, though.
b) I hope that, at some point, my library will be of such a quality
that it may be useful to others, and in that event I will release it.
Interoperability with other libraries is therefore a goal for me, and
a part of this is to let the user choose other vector/matrix types
than the ones provided by me.
Yes please. It would be great if you considered submitting it to Phobos.
I agree that a set of decent vector and matrix types should be put into
phobos. But the packages I mention above are perhaps more suited for a
dedicated scientific library, don't you think?
Right now I use D1, but I've been looking at the new phobos lately, and
I have to say that the combination D2+phobos2 is a very enticing one.
Writing array wrappers has never been easier than with the new "alias
this" feature. :)
If it weren't for the lack of a 64-bit compiler for Linux I would switch
immediately. As it is, I am seriously considering using the 32-bit one,
even though it just feels... wrong.
Actually, I tried compiling my lib with the latest 32-bit DMD2.
Immediately after pressing enter, memory usage went through the roof and
my computer became completely unresponsive. It took me forever to kill
the dmd process. I guess it has something to do with the heavy use of
templates. Has anyone else experienced this?
c) Often, for reasons of both efficiency and simplicity, it is
desirable to use arrays directly.
Yah.
My first question goes to those among you who do a lot of linear
algebra in D: Do you think supporting both library types and arrays
is worth the trouble? Or should I just go with one and be done with it?
If you go templated you don't need to explicitly support built-in arrays
- they'll just work.
Yeah, vectors work fine, especially now that D has built-in vector
operations. The problem is with matrices, as described in the following.
A user-defined matrix type would have opIndex(i,j) defined, and to
retrieve elements one would write m[i,j].
Yah.
However, the syntax for two-dimensional arrays is m[i][j], and this
means I have to put a lot of static ifs around my code, in order to
check the type every time I access a matrix.
No, that's not a two-dimensional array; it's an array of arrays. If you
want to make your lib work with arrays of arrays, you could easily build
a little wrapper arround it (e.g. JaggedMatrix).
OK, technically it may not be what is called a two-dimensional array.
But it's the closest we've got, no? Or perhaps that would be a
rectangular (static) array, but those still have the [][] indexing syntax.
This leads me to my second question, which is a suggestion for a
language change, so I expect a lot of resistance. :)
Would it be problematic to define m[i,j,...] to be equivalent to
m[i][j][...] for built-in arrays, so that arrays and user-defined
types could be used interchangeably?
(And, importantly, are there anyone but me who think they would
benefit from this?)
This wouldn't harm, but it would be a special case.
I know, and possibly an esoteric one, at that. But it would make it
easier to create drop-in array replacements.
-Lars
