Hi,
Firstly, sorry about the large amount of code in this question. OK, so I'm making a C++ wrapper for the gsl_matrix
struct and its associated functions. To make my approach a bit easier to extend to other kinds of gsl structs, I have a template base class called gsl_base
, which has a set of useful functions for getting at and checking the status of the underlying gsl struct when using the C++ wrapper class. So, gsl_base
looks like this:
namespace gsl{
template< typename T >
class gsl_base
{
public:
/// Get a pointer to the underlying GSL struct
///
/// Will not throw, can return NULL
inline T* ptr(){ return M_pGSLData; }
/// Get a const pointer to the underlying GSL struct
///
/// Will not throw, can return NULL
inline const T* const_ptr() const{ return M_pGSLData; }
/// Check if it's OK to try and manipulate the object
///
/// Will not throw
inline bool hasValue() const{ return M_pGSLData != 0; }
/// Check if it's not OK to try and manipulate the object
///
/// Will not throw
inline bool isNull() const{ return M_pGSLData == 0; }
protected:
gsl_base() : M_pGSLData( NULL ){}
gsl_base( T* original ) : M_pGSLData( original ){}
~gsl_base(){}
/// Set the value of the underlying GSL struct
///
/// Will not throw
inline void set_ptr( T* p ){ M_pGSLData = p; }
T* M_pGSLData;
};
}
And to make a simple C++ wrapper for gsl_matrix
I inherit from gsl_base
like this:
#pragma once
#include <stdexcept>
#include <cassert>
#include <iostream>
#include "../Common/number.h"
#include "../Common/gslstructbase.h"
#include <gsl/gsl_matrix.h>
namespace gsl
{
class realMatrix : public gsl::gsl_base< gsl_matrix > {
public:
typedef real value_type;
typedef size_t size_type;
typedef real* iterator;
typedef const real* const_iterator;
typedef real& reference;
typedef const real& const_reference;
typedef real* pointer;
typedef const real* const_pointer;
typedef ptrdiff_t difference_type;
struct M_matrix_element_type{
realMatrix::size_type row;
realMatrix::size_type col;
M_matrix_element_type( realMatrix::size_type r, realMatrix::size_type c ) :
row( r ), col( c ) {}
};
typedef M_matrix_element_type element_type;
/// Constructors
realMatrix();
realMatrix( size_type r, size_type c ) throw ( std::bad_alloc );
~realMatrix(){}
/// Print the addresses of the elements in the matrix (for debugging only!)
/// to be removed when problems with memory allocation are fixed
void printAddresses() const
{
std::cout << std::endl;
for ( size_t i = 0; i < M_uRows; ++i ){
for( size_t j = 0; j < M_uCols; ++j )
std::cout << M_pStart + i*M_uCols + j << " ";
std::cout << std::endl;
}
}
private:
pointer M_pStart;
pointer M_pFinish;
size_type M_uRows;
size_type M_uCols;
/// Set all private fields to correspond to a gsl_matrix structs values
/// m must not be NULL
///
/// Will not throw
inline void M_set_all_fields( const gsl_matrix* m )
{
assert( m != 0 );
M_pStart = m->data;
M_pFinish = M_pStart + m->size1 * m->size2;
M_uRows = m->size1;
M_uCols = m->size2;
}
};
}
real
is a typedef for double
defined in ../Common/number.h
. And the corresponding cpp file for reamMatrix
is quite simple and looks like:
#include "matrix.h"
#include "../Utils/mathutils.h"
////////////////////////////////////////////////////////////
namespace gsl{
////////////////////////////////////////////////////////////
realMatrix::realMatrix() : gsl_base(0), M_pStart(0), M_pFinish(0), M_uRows(0), M_uCols(0) {}
////////////////////////////////////////////////////////////
realMatrix::realMatrix( realMatrix::size_type r, realMatrix::size_type c )
throw( std::bad_alloc ) : gsl_base( gsl_matrix_alloc( atLeastOne( r ), atLeastOne( c ) ) )
{
if ( isNull() )
throw std::bad_alloc();
M_set_all_fields( this->const_ptr() );
}
////////////////////////////////////////////////////////////
} // End of gsl namespace
////////////////////////////////////////////////////////////
There are a bunch of other functions for doing various matrixy things, but these are the minimum needed to show the issue. OK, I also have an auto-test that was supposed to check a boolean equals that I had defined for gsl::realMatrix
, but that's all been stripped away now and all it does is instantiate two matrices and print the addresses of their elements. (for diagnosing why I was getting segfaults and crash dumps). So, my test looks like:
void MatrixTestSuite::MatricesAreEqual()
{
const int iRows = 10;
const int iCols = 5;
gsl::realMatrix m_1( iRows, iCols );
std::cout << "m_1 addresses = ";
m_1.printAddresses();
gsl::realMatrix mat_2( iRows, iCols );
std::cout << "m_1 addresses = ";
m_1.printAddresses();
std::cout << "m_2 addresses = ";
mat_2.printAddresses();
}
That's it. Instantiate a matrix, print its element addresses. Instantiate a second matrix, print the address of the first matrix's elements again and then the second's elements. The results are not what I expected. I get something like this:
m_1 addresses =
0x81bfd70 0x81bfd78 0x81bfd80 0x81bfd88 0x81bfd90
0x81bfd98 0x81bfda0 0x81bfda8 0x81bfdb0 0x81bfdb8
0x81bfdc0 0x81bfdc8 0x81bfdd0 0x81bfdd8 0x81bfde0
0x81bfde8 0x81bfdf0 0x81bfdf8 0x81bfe00 0x81bfe08
0x81bfe10 0x81bfe18 0x81bfe20 0x81bfe28 0x81bfe30
0x81bfe38 0x81bfe40 0x81bfe48 0x81bfe50 0x81bfe58
0x81bfe60 0x81bfe68 0x81bfe70 0x81bfe78 0x81bfe80
0x81bfe88 0x81bfe90 0x81bfe98 0x81bfea0 0x81bfea8
0x81bfeb0 0x81bfeb8 0x81bfec0 0x81bfec8 0x81bfed0
0x81bfed8 0x81bfee0 0x81bfee8 0x81bfef0 0x81bfef8
m_1 addresses =
0x1 0x9 0x11 0x19 0x21
0x29 0x31 0x39 0x41 0x49
0x51 0x59 0x61 0x69 0x71
0x79 0x81 0x89 0x91 0x99
0xa1 0xa9 0xb1 0xb9 0xc1
0xc9 0xd1 0xd9 0xe1 0xe9
0xf1 0xf9 0x101 0x109 0x111
0x119 0x121 0x129 0x131 0x139
0x141 0x149 0x151 0x159 0x161
0x169 0x171 0x179 0x181 0x189
m_2 addresses =
0x81bff08 0x81bff10 0x81bff18 0x81bff20 0x81bff28
0x81bff30 0x81bff38 0x81bff40 0x81bff48 0x81bff50
0x81bff58 0x81bff60 0x81bff68 0x81bff70 0x81bff78
0x81bff80 0x81bff88 0x81bff90 0x81bff98 0x81bffa0
0x81bffa8 0x81bffb0 0x81bffb8 0x81bffc0 0x81bffc8
0x81bffd0 0x81bffd8 0x81bffe0 0x81bffe8 0x81bfff0
0x81bfff8 0x81c0000 0x81c0008 0x81c0010 0x81c0018
0x81c0020 0x81c0028 0x81c0030 0x81c0038 0x81c0040
0x81c0048 0x81c0050 0x81c0058 0x81c0060 0x81c0068
0x81c0070 0x81c0078 0x81c0080 0x81c0088 0x81c0090
As you can see, the elements of the first matrix are invalidated by the instantiation of the second!?!? This is bugging the hell out of me, so any suggestions would be really gratefully received! Also, the command that's being used to compile this code is: g++ -c "gsl++/Source/gslpp/Matrix/matrix.cpp" -g -pedantic -Wall -o ../../Matrix_matrix.o -I. -I/usr/local/include -I/usr/include
Cheers.