gsl matrix

gsl matrix

[Juan Jose Gomez Cadenas]

Hello,
We have started to use gsl recently. Currently my group is involved in a
high energy physics experiment at CERN. We are developping the software
for the event reconstruction and display.

Our software is C++ and we use heavyly matrix operations. Thus I had the
idea to build a C++ wrapper class around the gsl matrices. It looks a good
idea to me (opinions?) since the adressing mechanism of the gsl matrices,
"a la valarray" should be fast and they come with a full BLAS interface.

I have been working in the first prototypes and found a couple of problems
which I would like to discuss.

The first one has to do with the absence in the gsl-0.9 library of a
function of the type: 

double * gsl_matrix_ptr(const gsl_matrix * m, const size_t i, const size_t
j);

The function (and all its cousins for the different data types) is
prototyped in the corresponding headers but not implemented. It must have
been supressed on purpose since the implementation is trivial. On the
other hand I would argue such a function is necessary for the "write
element" operator, which requires a pointer. We implement the operator
like this:

double& Rmatrix::operator () (int row, int col) const {
  double* elementPointer = gsl_matrix_ptr(gslMatrix, row-1, col-1);
  return *elementPointer;

}

Could this function be (re)-added to the package?

This is my provisional implementation of the function:

BASE *
FUNCTION (gsl_matrix, ptr) (const TYPE (gsl_matrix) * m,
			    const size_t i, const size_t j)
{

  if (gsl_check_range)
    {
      if (i >= m->size1)	/* size_t is unsigned, can't be negative
*/
	{
	  gsl_error ("first index out of range",
		     __FILE__, __LINE__, GSL_EINVAL) ;
	}
      else if (j >= m->size2)	/* size_t is unsigned, can't be negative
*/
	{
	  gsl_error ("second index out of range",
		     __FILE__, __LINE__, GSL_EINVAL) ;
	}
    }
    return (BASE *) (m->data + MULTIPLICITY * (i * m->tda + j));

}

Notice that I use the gsl_error function rather than the macros used in
similar functions (that would not work here).

The second "problem" that I have found appears when
compiling my wrapper class:

make[1]: Entering directory `/home/jj/Packages/PMatrix/v1/cmt'
Rebuilding ../i386_linux22/PMatrix_dependencies.make
make[1]: Leaving directory `/home/jj/Packages/PMatrix/v1/cmt'
make[1]: Entering directory `/home/jj/Packages/PMatrix/v1/cmt'
../i386_linux22/MatrixMachine.o
cd ../i386_linux22/; g++ -c -I"/home/jj/Packages/PMatrix/v1"
-I"../../../GSL/v1/src" -I"usr/local/include" -I"../../../NAGC/v6/src"
-I"/home/jj/Packages/NAGC/v6" -I"/home/jj/Software/NAGC/6.0/include"
-Df2cFortran -D__cplusplus -DHAVE_INLINE -fPIC -shared -D_GNU_SOURCE
-Dlinux -DLINUX -Dunix -pipe -ansi -Wall -W  -D_GNU_SOURCE -o
MatrixMachine.o                   ../src/MatrixMachine.cpp
In file included from
/usr/local/include/gsl/gsl_vector_complex_long_double.h:26,
                 from
/usr/local/include/gsl/gsl_matrix_complex_long_double.h:26,
                 from /usr/local/include/gsl/gsl_matrix.h:4,
                 from /home/jj/Packages/PMatrix/v1/PMatrix/PMatrix.h:11,
                 from ../src/MatrixMachine.cpp:10:
/usr/local/include/gsl/gsl_vector_long_double.h:58: warning: non-static
const member `const
gsl_vector_long_double {anonymous union}::vector' in class without a
constructorIn file included from
/usr/local/include/gsl/gsl_matrix_complex_long_double.h:26,
                 from /usr/local/include/gsl/gsl_matrix.h:4,
                 from /home/jj/Packages/PMatrix/v1/PMatrix/PMatrix.h:11,
                 from ../src/MatrixMachine.cpp:10:
/usr/local/include/gsl/gsl_vector_complex_long_double.h:60: warning:
non-static const member `const
gsl_vector_complex_long_double {anonymous union}::vector' in class without
a constructor
In file included from /usr/local/include/gsl/gsl_matrix.h:4,
                 from /home/jj/Packages/PMatrix/v1/PMatrix/PMatrix.h:11,
                 from ../src/MatrixMachine.cpp:10:
/usr/local/include/gsl/gsl_matrix_complex_long_double.h:59: warning:
non-static const member `const
gsl_matrix_complex_long_double {anonymous union}::matrix' in class without
a constructor
In file included from
/usr/local/include/gsl/gsl_vector_complex_double.h:26,
                 from
/usr/local/include/gsl/gsl_matrix_complex_double.h:26,
                 from /usr/local/include/gsl/gsl_matrix.h:5,
                 from /home/jj/Packages/PMatrix/v1/PMatrix/PMatrix.h:11,
                 from ../src/MatrixMachine.cpp:10:
/usr/local/include/gsl/gsl_vector_double.h:58: warning: non-static const
member `const gsl_vector
{anonymous union}::vector' in class without a constructor
In file included from
/usr/local/include/gsl/gsl_matrix_complex_double.h:26,
                 from /usr/local/include/gsl/gsl_matrix.h:5,
                 from /home/jj/Packages/PMatrix/v1/PMatrix/PMatrix.h:11,
                 from ../src/MatrixMachine.cpp:10:
/usr/local/include/gsl/gsl_vector_complex_double.h:60: warning: non-static
const member `const
gsl_vector_complex {anonymous union}::vector' in class without a
constructor
In file included from /usr/local/include/gsl/gsl_matrix.h:5,
                 from /home/jj/Packages/PMatrix/v1/PMatrix/PMatrix.h:11,
                 from ../src/MatrixMachine.cpp:10:
/usr/local/include/gsl/gsl_matrix_complex_double.h:59: warning: non-static
const member `const
gsl_matrix_complex {anonymous union}::matrix' in class without a
constructor
In file included from
/usr/local/include/gsl/gsl_vector_complex_float.h:26,
                 from
/usr/local/include/gsl/gsl_matrix_complex_float.h:26,
                 from /usr/local/include/gsl/gsl_matrix.h:6,
                 from /home/jj/Packages/PMatrix/v1/PMatrix/PMatrix.h:11,
                 from ../src/MatrixMachine.cpp:10:
/usr/local/include/gsl/gsl_vector_float.h:58: warning: non-static const
member `const gsl_vector_float
{anonymous union}::vector' in class without a constructor
In file included from
/usr/local/include/gsl/gsl_matrix_complex_float.h:26,
                 from /usr/local/include/gsl/gsl_matrix.h:6,
                 from /home/jj/Packages/PMatrix/v1/PMatrix/PMatrix.h:11,
                 from ../src/MatrixMachine.cpp:10:
/usr/local/include/gsl/gsl_vector

and so on. Is all the time the same type of warning. I have no idea how to
fix it.

Apologies for the lengthy mail and congratulations for the impressive work
done with gsl. Our idea is to incorporate as much as possible of the
library to our experiment's software.

best regards,
jj



Juan Jose Gomez Cadenas
Department of Atomic and Nuclear Physics and IFIC
University of Valencia

Edificio de Institutos de Paterna
Apto Correos 2085
46071 Valencia

telf. 96-3864500
fax-  96-3864583

e-mail             gomez@mail.cern.ch
                   gomez@hal.ific.uv.es

RE: gsl matrix

[J.J. Gomez Cadenas]

Dear Brian, all,
As suspected that cure of the problem was easy enough. Simply adding an
empty constructor within a #ifdef __cpluplus will do it.
I have know a version with this patch added to all the relevant header files
(there are quite a few!!!). I'd be pleased, of course, to send it to you, or
to committ it to the repository, or whatever. Please tell me what you
prefer.

Right now I have implemented a c++ real matrix class whith the basic
functionality (all the standar operations in the matrix ring). It works
fine, is cute, is fast and it's amazingly simple (just one private attribute
a pointer to a gsl_matrix!). Next steps will be to create a complex matrix
class (perhaps also integer) plus the corresponding classes for vectors and
some linear algebra. All this, of course would be also available for whoever
interested.

Best, jj

-----Mensaje original-----
De: Brian Gough [ mailto:bjg@network-theory.co.uk ]
Enviado el: viernes, 03 de agosto de 2001 19:44
Para: Juan Jose Gomez Cadenas
CC: gsl-discuss@sources.redhat.com
Asunto: Re: gsl matrix


Juan Jose Gomez Cadenas writes:
 > cd ../i386_linux22/; g++ -c -I"/home/jj/Packages/PMatrix/v1"
 > -I"../../../GSL/v1/src" -I"usr/local/include" -I"../../../NAGC/v6/src"
 > -I"/home/jj/Packages/NAGC/v6" -I"/home/jj/Software/NAGC/6.0/include"
 > -Df2cFortran -D__cplusplus -DHAVE_INLINE -fPIC -shared -D_GNU_SOURCE
 > -Dlinux -DLINUX -Dunix -pipe -ansi -Wall -W  -D_GNU_SOURCE -o
 > MatrixMachine.o                   ../src/MatrixMachine.cpp
 > /usr/local/include/gsl/gsl_vector_long_double.h:58: warning: non-static
 > const member `const
 > gsl_vector_long_double {anonymous union}::vector' in class without a
 > constructorIn file included from

I found that the error doesn't occur with g++-3.0, but does with g++ < 3.0.

I don't know if that means that the message was removed deliberately,
or if it is just missing from g++ 3.0 and will be put back in later.

Anyway I just thought I would mention it.

RE: gsl matrix

[Brian Gough]

J.J. Gomez Cadenas writes:
 > I am sending you two zipped files, which corresponds to the matrix and
 > vector directories. The fixes are in (lots of) the headers files. All what I
 > have done is to add code like this:

Thanks for the files.  I tried adding the constructors, in different
ways, but I always ended up with a warning message of some sort in
another place.  I think I am just going to leave it -- it is only a
problem with -W.  -Wall is ok on all the compilers I tried.  It was
worth trying anyway.

regards
Brian Gough

Re: gsl matrix

[Brian Gough]

Juan Jose Gomez Cadenas writes:
 > cd ../i386_linux22/; g++ -c -I"/home/jj/Packages/PMatrix/v1"
 > -I"../../../GSL/v1/src" -I"usr/local/include" -I"../../../NAGC/v6/src"
 > -I"/home/jj/Packages/NAGC/v6" -I"/home/jj/Software/NAGC/6.0/include"
 > -Df2cFortran -D__cplusplus -DHAVE_INLINE -fPIC -shared -D_GNU_SOURCE
 > -Dlinux -DLINUX -Dunix -pipe -ansi -Wall -W  -D_GNU_SOURCE -o
 > MatrixMachine.o                   ../src/MatrixMachine.cpp
 > /usr/local/include/gsl/gsl_vector_long_double.h:58: warning: non-static
 > const member `const
 > gsl_vector_long_double {anonymous union}::vector' in class without a
 > constructorIn file included from

I found that the error doesn't occur with g++-3.0, but does with g++ < 3.0.

I don't know if that means that the message was removed deliberately,
or if it is just missing from g++ 3.0 and will be put back in later.

Anyway I just thought I would mention it.

RE: gsl matrix

[Brian Gough]

J.J. Gomez Cadenas writes:
 > As suspected that cure of the problem was easy enough. Simply adding an
 > empty constructor within a #ifdef __cpluplus will do it.
 > I have know a version with this patch added to all the relevant header files
 > (there are quite a few!!!). I'd be pleased, of course, to send it to you, or
 > to committ it to the repository, or whatever. Please tell me what you
 > prefer.

Ok, if you send the patches to me I will add them.

Re: gsl matrix

[Brian Gough]

Juan Jose Gomez Cadenas writes:
 > The first one has to do with the absence in the gsl-0.9 library of a
 > function of the type: 
 > 
 > double * gsl_matrix_ptr(const gsl_matrix * m, const size_t i, const size_t
 > j);
 >
 > Could this function be (re)-added to the package?

Good point, I've now added them in the CVS version.  I haven't tested
them so you will need to check that they work correctly.

 > Notice that I use the gsl_error function rather than the macros used in
 > similar functions (that would not work here).

I used the GSL_ERROR_NULL macro to return a NULL pointer for
out-of-range arguments.  I think this should work ok.

 > The second "problem" that I have found appears when
 > compiling my wrapper class:
 > 
 > /usr/local/include/gsl/gsl_vector_long_double.h:58: warning: non-static
 > const member `const
 > gsl_vector_long_double {anonymous union}::vector' in class without a
 > constructor

Maybe we'll need to have a empty constructor inside an #ifdef
__cplusplus .. #endif to get rid of the warning.