newbie dynamic memory problem.

[laclac01]

I wanted to play around and see if I could make a simple little set of
routines that will do simple complex matrix math. I am not very far in
to but already it is crashing. I am sure I am not doing something
correctly so I was hoping someone smarter than I (meaning almost
everyone) could take a look at what I have and tell me what's wrong.

struct complexMatrix
{
float r;
float i;
};

#include<iostream.h>
#include "types.h"
void CmxAlloc(complexMatrix** mat,int rows, int cols);
void CmxFree(complexMatrix** mat);

int main()
{
complexMatrix** my;

CmxAlloc(my,10,10);

my[0][0].i=6; //<- this causes the program to crash




return 0;
}
//*******************************************************
void CmxAlloc(complexMatrix** mat, int rows, int cols)
{
mat = new complexMatrix*[rows];
for(int i = 0; i < rows; i++)
mat = new complexMatrix[cols];

}
//*******************************************************
void CmxFree(complexMatrix** mat)
{
delete [] mat;
}

 

[Victor Bazarov]

I wanted to play around and see if I could make a simple little set of
routines that will do simple complex matrix math. I am not very far in
to but already it is crashing. I am sure I am not doing something
correctly so I was hoping someone smarter than I (meaning almost
everyone) could take a look at what I have and tell me what's wrong.

struct complexMatrix
{
float r;
float i;
};

#include<iostream.h>
#include "types.h"
void CmxAlloc(complexMatrix** mat,int rows, int cols);
void CmxFree(complexMatrix** mat);

Adjust the prototypes accordingly. See below.

int main()
{
complexMatrix** my;

CmxAlloc(my,10,10);

my[0][0].i=6; //<- this causes the program to crash




return 0;
}
//*******************************************************
void CmxAlloc(complexMatrix** mat, int rows, int cols)

You're passing the pointer 'mat' by value. You change it inside this
function. Nobody outside the function knows anything about the changes.
Pass by reference:

void CmxAlloc(complexMatrix**& mat, int rows, int cols)

Now, whatever you assign to 'mat' will be visible to the caller of this
function.

{
mat = new complexMatrix*[rows];
for(int i = 0; i < rows; i++)
mat = new complexMatrix[cols];

}

Of course it would be better to make your function like this:

complexMatrix ** CmxAlloc(int rows, int cols)
{
complexMatrix ** mat = new complexMatrix*[rows];
for (int i = 0; i < rows; i++)
mat = new complexMatrix[cols];
return mat;
}

then there would be no need to pass by reference.

//*******************************************************
void CmxFree(complexMatrix** mat)
{
delete [] mat;

You're only deleting the array of pointers. Each pointer in that array
remains unfreed -- you leak memory. Need to pass the 'rows' and do the
inverse of 'CmxAlloc':

for (int i = rows - 1; i >= 0; --i)
delete[] mat;
delete[] mat;

}

V

 

[Esteban Osses]

your function CmxAlloc takes the complexMatrix** parameter by value, so
you allocated the space for your complexMatrix, but you never assigned
the direction of that matrix to `my'.

 

[Kai-Uwe Bux]

I wanted to play around and see if I could make a simple little set of
routines that will do simple complex matrix math. I am not very far in
to but already it is crashing. I am sure I am not doing something
correctly so I was hoping someone smarter than I (meaning almost
everyone) could take a look at what I have and tell me what's wrong.

struct complexMatrix
{
float r;
float i;
};

#include<iostream.h>
#include "types.h"
void CmxAlloc(complexMatrix** mat,int rows, int cols);
void CmxFree(complexMatrix** mat);

int main()
{
complexMatrix** my;

CmxAlloc(my,10,10);

my[0][0].i=6; //<- this causes the program to crash




return 0;
}
//*******************************************************
void CmxAlloc(complexMatrix** mat, int rows, int cols)

make that:

void CmxAlloc(complexMatrix**& mat, int rows, int cols)

{
mat = new complexMatrix*[rows];
for(int i = 0; i < rows; i++)
mat = new complexMatrix[cols];

}
//*******************************************************
void CmxFree(complexMatrix** mat)
{
delete [] mat;

That leaks big time: none of your rows will actually be deallocated.

}

BTW, do you really think managing memory explicitly yourself is the right
way to go?


Best

Kai-Uwe Bux

 

[Frank Chang]

Kai-Uwe Bux, I agree with you that "managing memory explicity yourself
is not the right way to go in this case." I have a question of style,
couldn't you write CmxAlloc as either

void CmxAlloc(complexMatrix**& mat, int rows, int cols)

or

void CmxAlloc(complexMatrix*** mat, int rows, int cols) ?

If so, which is the preferred way? Thank you.

 

[Victor Bazarov]

Kai-Uwe Bux, I agree with you that "managing memory explicity yourself
is not the right way to go in this case." I have a question of style,
couldn't you write CmxAlloc as either

void CmxAlloc(complexMatrix**& mat, int rows, int cols)

or

void CmxAlloc(complexMatrix*** mat, int rows, int cols) ?

If so, which is the preferred way? Thank you.

Sorry to barge in, my preferred way would be

template<class T> T* Alloc(int howmany)
{
return new T[howmany];
}

template<class T> T** Alloc(int rows, int cols)
{
T **tt = Alloc<T*>(rows);
while (rows)
tt[--rows] = Alloc<T>(cols);
return tt;
}

template<class T> void Free(T* tt)
{
delete[] tt;
}

template<class T> void Free(T** tt, int rows)
{
for (int i = 0; tt && i < rows; ++i)
Free(tt);
Free(tt);
}
...
complexMatrix **mat = Alloc<complexMatrix>(rows, cols);
...
Free(mat, rows);

Feel free to poke holes in it.

V

 

[Frank Chang]

Victor, I haven't tested your code on a large symmetric Complex matrix
There are no holes in your code , except why don't you use a template
class:

template <class Type>
class Matrix {
public:
Matrix(int rows, int cols);
Matrix(const Matrix& iM);
virtual ~Matrix();
Matrix& operator=(const Matrix& iM);

protected:
Type** Alloc(int rows, int cols);
Type** Free(T** tt, int rows)

int nrows;
int ncols;
Type** mat
}

 

[Victor Bazarov]

Victor, I haven't tested your code on a large symmetric Complex matrix
There are no holes in your code , except why don't you use a template
class:

template <class Type>
class Matrix {
public:
Matrix(int rows, int cols);
Matrix(const Matrix& iM);
virtual ~Matrix();
Matrix& operator=(const Matrix& iM);

protected:
Type** Alloc(int rows, int cols);
Type** Free(T** tt, int rows)

int nrows;
int ncols;
Type** mat
}

The code was posted by the original poster. I didn't propose, nor do I
endorse, using a two-dimensional array for a matrix of complex numbers.
It's his/her choice, I am just giving possible corrections for the simple
task at hand: allocating two-dim arrays using 'new[]' in a function.

Why don't you ask the original poster about using class Matrix?

V

 

[Frank Chang]

Victor, Thank you for clarifying the design issue. I have tried many
times asking questions to the original poster about his/her symmetric
matrix and memory constraints and he/she never responds. There are so
many interesting questions popping up all the time in this newsgroup it
is impossible to keep track of all the questions.

 

[Kai-Uwe Bux]

What is the best way to do this??? I am not sure I am familiar with the
matrix class. What isn't a 2d array the best way to do this?

A 2d array, i.e., something like

double** matr

is *bad*: the object does not know anything about its size. Thus you have to
keep track of its number of rows and columns separately so that you clean
up the mess when you want to get rid of the matrix. Those nitty gritty
details can be really hard to get right.

This is actually just a special case of arrays being bad. Compare double**
to:

typedef std::vector< double > MathVector;
typedef std::vector< MathVector > MathMatrix;


These objects will manage the allocated memory automatically for you.

Also, these objects know about their size. Thus, you could write things like

MathVector operator+ ( MathVector const & v, MathVector const & u ) {
if ( v.size() != u.size() ) {
throw; // should throw something specific;
}
MathVector result ( v.size(), 0.0 );
for ( unsigned long i = 0; i < v.size(); ++i ) {
result = v + u;
}
return( result );
}

or

std:stream & operator<< ( std:stream & str,
MathVector v ) {
str << '(';
for ( unsigned int i = 0; i < v.size(); ++i ) {
str << ' ' << v;
}
return( str << " )" );
}



This is just meant as a hint. I would not really recommend using MathMatrix
as a representation of a matrix since it does not guarantee that all
Vectors have the same size. Thus, one should really wrap this into a class
that enforces rectangular shape.

Finally, since we were talking about matrix classes, a word of warning:

a) there are powerful linear algebra libraries out there whose design and
performance are very hard to beat.

b) numerical analysis is *hard*: numbers (float or double) on a computer are
very different from numbers in your linear algebra text book. Thus, what
looks like a correct algorithm might turn out to be completely misguided,
e.g., because of rounding errors. For instance, determinants close to 0 are
usually a problem. These things are best left to experts; and this means
experts in numerical analysis not in programming.


Thus: roll your own code only for educational purposes.


Best

Kai-Uwe Bux

 

[Kai-Uwe Bux]

I wanted to play around and see if I could make a simple little set of
routines that will do simple complex matrix math. I am not very far in
to but already it is crashing. I am sure I am not doing something
correctly so I was hoping someone smarter than I (meaning almost
everyone) could take a look at what I have and tell me what's wrong.

struct complexMatrix
{
float r;
float i;
};

In addition to other remarks somewhere in this thread, I would like to
suggest you have a look into

std::complex<>


Best

Kai-Uwe Bux

 

[Frank Chang]

Victor, Yes , "there are many implementations of matrix on the 'Net. I
think lacla... should look at the Boost matrix class. The Boost
distribution is easy to build, install and use(less than 3 hours). Many
of the matrix class libraries on the market and 'Net are overly complex
and hard to use because they require that somebody have a thorough
understanding of linear algebra in order to properly uses these matrix
classes.

Victor,you state "A 2d array wastes a bit more memory than it
should.". This is only true for certain classes of matrices . In many
applications, all the elements of a 2-d array are used.

Also, you state "There is no automatic disposal of the memory when
you're done with it." In many C++ applications , the onus for the
disposal of the memory is on the programmer who needs to write the
correct destructor. There is no automatic panacea for preventing memory
leaks. Unmanaged C++ does not use garbage collection unlike other
languages such as managed C++, C#, Java, Scheme..

 

[Frank Chang]

Kai-Uwe Bux, Yes, IN STL , you can use a vector of vectors to
represent a matrix. As you state, the STL vectors can take the
responsibility for the dynamic memory management. However, the
2d array, double** matr, while challenging to code , may allow random
access to matrix(i, j) to occur faster than with a vector of vectors.
This may make a difference in a very large matrix computation.

 

[Kai-Uwe Bux]

Kai-Uwe Bux, Yes, IN STL , you can use a vector of vectors to
represent a matrix. As you state, the STL vectors can take the
responsibility for the dynamic memory management. However, the
2d array, double** matr, while challenging to code , may allow random
access to matrix(i, j) to occur faster than with a vector of vectors.
This may make a difference in a very large matrix computation.

Generally, the notion that std::vector< T > has a performance penalty
compared to T* is, in my opinion, misguided. However, this is not your
point, as you are specifically thinking about chaining calls to operator[].

Nonetheless, I doubt that std::vector< std::vector< T > > has a disadvantage
compared to T**. In its most simple incarnation, an implementation of
std::vector<T> would, maybe, start like so:

template < typename T >
class vector {
private:

T * data;
std::site_t length;

public:

...

T& operator[] ( std::size_t pos ) {
return( data[pos] );
}

...

}; // class vector<>

Note that the operator[] call will most definitely be inlined directly by
the compiler. Also note that the data-field starts the struct. That means:
if the compiler is sufficiently simple-minded about memory layout, a
pointer to a std::vector<T> object will point directly to the T* field
inside the object. Thus chaining operator[]-calls from
std::vector< std::vector<T> >
incurs no overhead as compared to working with T**. Very likely, the
compiler will actually produce identical code.


Also, I am *not* advocating std::vector< std::vector< double > > as an
approach top implement a high-performance linear algebra library.



Best

Kai-Uwe Bux

 

[Ram]

Kai-Uwe Bux, Yes, IN STL , you can use a vector of vectors to

represent a matrix. As you state, the STL vectors can take the
responsibility for the dynamic memory management. However, the
2d array, double** matr, while challenging to code , may allow random
access to matrix(i, j) to occur faster than with a vector of vectors.
This may make a difference in a very large matrix computation.

I agree with Kai-Uwe Bux that there may not be performance penalty
while accessing an element from a chained vector. However for another
reason, chaining of vectors can be an overkill in situations like
Matrix where the number of columns in each row is same. In this case,
each vector implementing a row also keeps an extra information which is
its size (no. of columns). This can lead to wastage of a large amount
of memory.

Ram

 

[Frank Chang]

Kai-Uwe Bux, I agree wholeheartedly with everything you say, but I
looked in the vector template class in VC7 and found the following
code:

reference operator[](size_type _Off)
{ // subscript mutable sequence
return (*(begin() + _Off));
}

where begin() is defined as

iterator begin()
( // return iterator for beginning of mutable sequence
return (iterator(_ITER_BASE(_Myvec.begin())));
}

So, in the MSVC 7.1 STL implementation, chaining operator[] operators
may incur some overhead as compared to using T**. In any case, you have
already stated that you are not advocating std::vector< std::vector
<double>> as an approach to implement a high-performance linear algebra
library. I am experimenting with the most recent BOOST distribution
now. What is your opinion about the Boost matrix template C++ class?
Thank you for taking the time on a weekend to discuss this topic.

 

[Victor Bazarov]

[..]
Victor,you state "A 2d array wastes a bit more memory than it
should.". This is only true for certain classes of matrices . In many
applications, all the elements of a 2-d array are used.

If I were to implement a matrix (assuming it's fully represented, not
a band matrix or a profile-based scheme), I'd still allocate a plain
array to represent the elements. It's better than allocate N pointers
and then make each of them point to M data elements. Think about it.
Only _one_ allocation of N*M elements or N+M allocations? I am talking
about the overhead dynamic memory allocation carries.

Also, you state "There is no automatic disposal of the memory when
you're done with it." In many C++ applications , the onus for the
disposal of the memory is on the programmer who needs to write the
correct destructor.

If you don't give your 2d array to be handled in some class, then I
ask you: what destructor? There is no destructor here. OTOH, if the
memory is the full responsibility of a 'matrix' class, then yes, the
destructor of a 'matrix' object will dispose of the memory.

The OP's implementation was a raw 2d array of some structs. I was
just pointing out that wrapped in a class it would be easier to deal
with.

[..]

V

 

[Frank Chang]

laclac01:

What is the best way to do this??? I am not sure I am familiar with the
matrix class. What isn't a 2d array the best way to do this?

I think Victor Bazorov is saying that a 2d matrix is okay as long as
you create a template C++ 2D matrix class with all the canonical member
functions (constructor, copy constructor, destructor, assignment
operator) instead of a struct. It appears that between you and Victor,
the basic pieces of this class have already been written. If you wish,
post your 2D matrix class on the bulletin board after you have written
it and I am sure people will give you helpful comments. Alternatively,
since STL does not have a matrix class, why don't you download the
latest Boost(www.boost.org) distribution and look at Boost's matrix
class(matrix.hpp)?

 

[Kai-Uwe Bux]

Kai-Uwe Bux, I agree wholeheartedly with everything you say, but I
looked in the vector template class in VC7 and found the following
code:

reference operator[](size_type _Off)
{ // subscript mutable sequence
return (*(begin() + _Off));
}

where begin() is defined as

iterator begin()
( // return iterator for beginning of mutable sequence
return (iterator(_ITER_BASE(_Myvec.begin())));
}

So, in the MSVC 7.1 STL implementation, chaining operator[] operators
may incur some overhead as compared to using T**.

Do not judge from the source: sometimes those macros involve compiler magic.
Keep in mind this is the library from your vendor. Very likely it is
finetuned to work exceptionally well in conjunction with the optimization
routines of the compiler. Sometimes I try, just for fun, to beat the
performance of the library that ships with g++, but I rarely ever succeed.


Why not just test it:

// array vs vector
// ===============

#include <cstdlib>
#include <ctime>
#include <vector>
#include <iostream>

unsigned long runs = 10000;
unsigned long rows = 314;
unsigned long cols = 218;

typedef unsigned long number;
typedef number * RawVector;
typedef RawVector * RawMatrix;
typedef std::vector< number > NumVector;
typedef std::vector< NumVector > NumMatrix;

number run_array ( void ) {
RawMatrix matrix = new RawVector [ rows ];
for ( unsigned long i = 0; i < rows; ++i ) {
matrix = new number [ cols ];
}

for ( unsigned long r = 0; r < rows; ++r ) {
for ( unsigned long c = 0; c < cols; ++c ) {
matrix[r][c] = r*c;
}
}

std::clock_t ticks = std::clock();
number sum = 0;
for ( unsigned int i = 0; i< runs; ++i ) {
for ( unsigned long r = 0; r < rows; ++r ) {
for ( unsigned long c = 0; c < cols; ++c ) {
sum += matrix[r][c];
}
}
}
ticks = std::clock() - ticks;
std::cout << "array: " << ticks << '\n';

for ( unsigned long i = 0; i < rows; ++i ) {
delete [] matrix;
}
delete [] matrix;

// we return the result of the computation so that the
// loop cannot be optimzed away.
return ( sum );
}

number run_vector ( void ) {
NumMatrix matrix ( rows, NumVector( cols, 0 ) );

for ( unsigned long r = 0; r < rows; ++r ) {
for ( unsigned long c = 0; c < cols; ++c ) {
matrix[r][c] = r*c;
}
}

std::clock_t ticks = std::clock();
number sum = 0;
for ( unsigned int i = 0; i< runs; ++i ) {
for ( unsigned long r = 0; r < rows; ++r ) {
for ( unsigned long c = 0; c < cols; ++c ) {
sum += matrix[r][c];
}
}
}
ticks = std::clock() - ticks;
std::cout << "vector: " << ticks << '\n';

// we return the result of the computation so that the
// loop cannot be optimzed away.
return ( sum );
}


int main ( void ) {
number sum = 0;
std::srand( 3458612 );
while ( true ) {
if ( std::rand() % 2 == 1 ) {
sum = run_vector();
} else {
sum = run_array();
}
}
}


Be sure to compile with optimization turned on. I got:
news_group> a.out
array: 1350000
array: 1410000
array: 1400000
vector: 1380000
vector: 1370000
vector: 1430000
array: 1460000
array: 1460000
vector: 1410000
array: 1420000
array: 1460000
vector: 1420000

In any case, you have
already stated that you are not advocating std::vector< std::vector
<double>> as an approach to implement a high-performance linear algebra
library. I am experimenting with the most recent BOOST distribution
now. What is your opinion about the Boost matrix template C++ class?

Never tried, no opinion.


Best

Kai-Uwe Bux

 

[Frank Chang]

Kai-Uwe Bux, I downloaded a recent GNU g++ STL release and here is the
vector source for your STL chaining implementation of the 2D matrix:

/**
* @brief Subscript access to the data contained in the vector.
* @param n The element for which data should be accessed.
* @return Read/write reference to data.
*
* This operator allows for easy, array-style, data access.
* Note that data access with this operator is unchecked and
out_of_range
* lookups are not defined. (For checked lookups see at().)
*/
reference operator[](size_type __n) { return *(begin() + __n); }

/**
* Returns a read/write iterator that points to the first element in
the
* vector. Iteration is done in ordinary element order.
*/
iterator begin() { return iterator (_M_start); }

/**
* @brief Provides access to the data contained in the vector.
* @param n The element for which data should be accessed.
* @return Read/write reference to data.
*
* This function provides for safer data access. The parameter is
first
* checked that it is in the range of the vector. The function
throws
* out_of_range if the check fails.
*/
reference at(size_type __n)
{ _M_range_check(__n); return (*this)[__n]; }

void _M_range_check(size_type __n) const {
if (__n >= this->size())
__throw_out_of_range("vector");
}

As you can see, this is identical to the MSVC7.1 STL code I posted the
other day. So, your argumemt about STL compiler differences doesn't
apply to this case. Yes, I agree each vendor uses compiler
optimizations. The MSVC7.1 STL header files list all the compiler
optimizations. As for your tests, I haven't ran them yet. At first
glance, the matrices in your test are not large and it appears you are
not using the standard linear algebra / matrix benchmarks which would
be a better test for random access of 2D matrices.