I am programming a Matrix class to be used like a primitive type (overloaded operators, identity, transpose, the whole shebang) and one of the constructors _may_ cause an Access Violation (or Segmentation Fault, if you prefer the old Linux SIGSEG) if the class's user is foolish. Is there a way to catch this errors?

Secondly, would this class be of any use to anyone? I'll gladly distribute it all around (OpenSourced, of course). It is 'templatized'. I'm unsure if Phobos or Tango would want such class if I submited it.

Silverling wrote:
> I am programming a Matrix class to be used like a primitive type (overloaded operators, identity, transpose, the whole shebang) and one of the constructors _may_ cause an Access Violation (or Segmentation Fault, if you prefer the old Linux SIGSEG) if the class's user is foolish. Is there a way to catch this errors?
> 
> Secondly, would this class be of any use to anyone? I'll gladly distribute it all around (OpenSourced, of course). It is 'templatized'. I'm unsure if Phobos or Tango would want such class if I submited it.

What are the specifics of your class?  Is it an arbitrary MxN sized matrix?  Will it use vendor-provided accelerated BLAS libraries where available?  Is it parameterized on element type?  Storage format? Can it handle different storage schemes (e.g. sparse formats like CSC, CSR, banded, symmetric).  If it's got all that, then I'm definitely interested.  A port of something like Boost::ublas would be great. (http://www.boost.org/libs/numeric/ublas/doc/index.htm).

I've started on a multidimensional array class that can be used as a matrix if desired.  It's based loosely on NumPy's multidimensional arrays.  http://www.dsource.org/projects/multiarray.html .   It's hooked up to BLAS and LAPACK for some basic operations like multiplication, linear system solving, and least squares problems.   Doesn't do sparse at all though.

Oskar Linde also posted a multidim array struct a while back.  Same basic philosophy, but using a struct instead of a class.  I'm probably going to switch over to using a struct eventually too, once D gets const ref.

--bb

> Is it an arbitrary MxN sized matrix?
Yes. It is resizeable after creation for purposes of multiplication and transpose.
> Will it use vendor-provided accelerated BLAS libraries where available?
Never heard of BLAS, but my I don't see why one couldn't change the code to use such libraries. Anyway, I'm keeping my module's dependencies down (currently depends only on std.string).

> Is it parameterized on element type?
I'm unsure of what you mean, but as I said before it is a template class. The matrices values could have any type, even other classes.

> Storage format?
Currently Type[row][col].

> Can it handle different storage schemes (e.g. sparse formats like CSC, CSR, banded, symmetric).
Not yet, but it can be altered to support them, probably not by me. I'll use it primarily to calculate view frustums (which usually don't have a lot of '0' justifying the current storage scheme) and vectorial calculus.
Implementing such schemes would only need to override the opIndex and opIndexAssign. The current implementation still accesses the matrix's data directly, but that can be easily changed (which I will, due to your idea).

I've started this module because I wanted to learn operator overloading on D. I'm having a small issue on overloading the opMul. I templatized it but I need a specialization to multipliy by a matrix. DMD is not accepting
void opMul(T:Matrix)(T mult)
meanwhile I've shifted my attention towards inverse matrix, determinant and equation system solving.
Later on I'll use it for Hamming coding (to port a script of mine from MatLab), where I'll test the matrix's behavior with type 'bool', although for memory purposes, it'll use bit arrays or (a lot of) masking.

I'll send the source if you'd like to take a peek.

Reply to silverling,

> I am programming a Matrix class to be used like a primitive type
> (overloaded operators, identity, transpose, the whole shebang) and one
> of the constructors _may_ cause an Access Violation (or Segmentation
> Fault, if you prefer the old Linux SIGSEG) if the class's user is
> foolish. Is there a way to catch this errors?

Put this in somewhere at global scope and all seg-v's will throw an Error
by passing signal a different function you can make it do other things.

	// stubs for access to unix signal handeling alias void function(int) sighandler_t;

extern (C) sighandler_t signal(int signum, sighandler_t handler);

const int SIGHUP	= 1;	// Term Hangup detected on controlling terminal or death of controlling process
const int SIGINT	= 2;	// Interrupt from keyboard
const int SIGQUIT	= 3;	// Quit from keyboard
const int SIGILL	= 4;	// Illegal Instruction
const int SIGABRT	= 6;	// Abort signal from abort(3)
const int SIGFPE	= 8;	// Floating point exception
const int SIGKILL	= 9;	// Kill signal
const int SIGSEGV	= 11;	// Invalid memory reference
const int SIGPIPE	= 13;	// Broken pipe: write to pipe with no readers
const int SIGALRM	= 14;	// Timer signal from alarm(2)
const int SIGTERM	= 15;	// Termination signal

static this()
{
		// trap seg-v and report as error (not sure this is working)
	signal(SIGSEGV,function void(int i){throw new Error("SEGV");});
}

BCS Wrote:

> Reply to silverling,
> 
> > I am programming a Matrix class to be used like a primitive type (overloaded operators, identity, transpose, the whole shebang) and one of the constructors _may_ cause an Access Violation (or Segmentation Fault, if you prefer the old Linux SIGSEG) if the class's user is foolish. Is there a way to catch this errors?
> 
> Put this in somewhere at global scope and all seg-v's will throw an Error by passing signal a different function you can make it do other things.
> 
> 	// stubs for access to unix signal handeling
> alias void function(int) sighandler_t;
> 
> extern (C) sighandler_t signal(int signum, sighandler_t handler);
> 
> const int SIGHUP	= 1;	// Term Hangup detected on controlling terminal or
> death of controlling process
> const int SIGINT	= 2;	// Interrupt from keyboard
> const int SIGQUIT	= 3;	// Quit from keyboard
> const int SIGILL	= 4;	// Illegal Instruction
> const int SIGABRT	= 6;	// Abort signal from abort(3)
> const int SIGFPE	= 8;	// Floating point exception
> const int SIGKILL	= 9;	// Kill signal
> const int SIGSEGV	= 11;	// Invalid memory reference
> const int SIGPIPE	= 13;	// Broken pipe: write to pipe with no readers
> const int SIGALRM	= 14;	// Timer signal from alarm(2)
> const int SIGTERM	= 15;	// Termination signal
> 
> static this()
> {
> 		// trap seg-v and report as error (not sure this is working)
> 	signal(SIGSEGV,function void(int i){throw new Error("SEGV");});
> }
> 
> 


This is great! Thanks! Now I only need to throw an exception. Unfortunatly, it will only work on linux :( meaning the module is not portable.

Reply to silverling,

> BCS Wrote:
> 
>> Reply to silverling,
>> 
>>> I am programming a Matrix class to be used like a primitive type
>>> (overloaded operators, identity, transpose, the whole shebang) and
>>> one of the constructors _may_ cause an Access Violation (or
>>> Segmentation Fault, if you prefer the old Linux SIGSEG) if the
>>> class's user is foolish. Is there a way to catch this errors?
>>> 
> This is great! Thanks! Now I only need to throw an exception.
> Unfortunatly, it will only work on linux :( meaning the module is not
> portable.
> 

so just tell the win32 user to not be foolish <g>

Silverling wrote:
>> Is it an arbitrary MxN sized matrix?
> Yes. It is resizeable after creation for purposes of multiplication and transpose.

>> Will it use vendor-provided accelerated BLAS libraries where available?
> Never heard of BLAS, but my I don't see why one couldn't change the code to use such libraries. Anyway, I'm keeping my module's dependencies down (currently depends only on std.string).

BLAS is good if you're going to be working with big matrices.  But for a 4x4 view mat, there's not much point.

>> Is it parameterized on element type?
> I'm unsure of what you mean, but as I said before it is a template class. The matrices values could have any type, even other classes.

Yeh, that's what I mean.  "parameterized on type" == The type of the elements is a template parameter.

>> Storage format?
> Currently Type[row][col].

Got it.  You may be better off with Type[row*col].  Type[row][col] is an array of col pointers to 1D arrays of rows, rather than densely packed memory.

> 
>> Can it handle different storage schemes (e.g. sparse formats like CSC, CSR, banded, symmetric).
> Not yet, but it can be altered to support them, probably not by me. I'll use it primarily to calculate view frustums (which usually don't have a lot of '0' justifying the current storage scheme)

For 3D linear algebra (or 4D homogeneous linalg) there's helix. http://www.dsource.org/projects/helix  It's pretty decent -- even includes a polar decomposition routine.  It didn't have 2D classes, which I missed, so I added those to my own copy.
There've been a few other folks to implement 3D linalg stuff in D as well.  Don't have any links for those though.  Anyway, for 3D stuff I think using a generic MxN matrix class is overkill.

> and vectorial calculus.
> Implementing such schemes would only need to override the opIndex and opIndexAssign. 
> The current implementation still accesses the matrix's data directly,
> but that can be easily changed (which I will, due to your idea).

But if you have different storage schemes you want to be able to do things like sparse_mat * dense_vec,  or dense_mat + upper_triangular_mat.   It's a multiple dispatch problem so it requires some thinking.  The approach I've seen is to

> I've started this module because I wanted to learn operator overloading on D. I'm having a small issue on overloading the opMul. I templatized it but I need a specialization to multipliy by a matrix. DMD is not accepting
> void opMul(T:Matrix)(T mult)

My understanding is that if you use specialization in D it disables IFTI (the thing that allows to to call a template function without specifying parameters).  So to use that you'd need to call it explicitly as A.opMul!(typeof(B))(B).  Not so nice.  The workaround is generally to use some static ifs or static asserts inside the template instead.
    void opMul(T)(T mult) {
         static assert(is(T:Matrix), "Mult only works for things derived from Matrix");
    }


--bb

Bill Baxter wrote:

> But if you have different storage schemes you want to be able to do things like sparse_mat * dense_vec,  or dense_mat + upper_triangular_mat.   It's a multiple dispatch problem so it requires some thinking.  The approach I've seen is to

Oops forgot to compelete that thought.  Take a look at things like Boost::ublas or MTL (Matrix Template Library) if you want some ideas how to deal with different storage mechanisms.

--bb

Bill Baxter wrote:
> Silverling wrote:
>> Bill Baxter wrote:
>>> Storage format?
>> Currently Type[row][col].
> 
> Got it.  You may be better off with Type[row*col].  Type[row][col] is an array of col pointers to 1D arrays of rows, rather than densely packed memory.

No, it's definitely densely packed memory. Static arrays always put their elements directly where you put the array itself, "inline".

Only dynamic (Type[]) and associative (Type[Type2]) arrays use "hidden" pointers.

Frits van Bommel wrote:
> Bill Baxter wrote:
>> Silverling wrote:
>  >> Bill Baxter wrote:
>>>> Storage format?
>>> Currently Type[row][col].
>>
>> Got it.  You may be better off with Type[row*col].  Type[row][col] is an array of col pointers to 1D arrays of rows, rather than densely packed memory.
> 
> No, it's definitely densely packed memory. Static arrays always put their elements directly where you put the array itself, "inline".
> 
> Only dynamic (Type[]) and associative (Type[Type2]) arrays use "hidden" pointers.

Ok, but he said it was going to be resizeable.  Does that change your answer?

--bb