difficulty understanding the data structure in meschach library

Question

I am reading the source code of the MESCHACH library for matrix and vector computations. Here is the data structure used in the library:

/* matrix definition */
typedef struct  {
        unsigned int    m, n;
        unsigned int    max_m, max_n, max_size;
        Real    **me,*base; /* base is base of alloc'd mem */
        } MAT;

In here, what is the use of *base in the structure ? I mean **me is the one containing the matrix and the upper values contain its dimensions, what does the base hold.

Code for allocating memory for the matrix:

MAT *m_get(int m, int n)
{
   MAT  *matrix;
   int  i;
   
   if (m < 0 || n < 0)
     error(E_NEG,"m_get");

   if ((matrix=NEW(MAT)) == (MAT *)NULL )
     error(E_MEM,"m_get");
   else if (mem_info_is_on()) {
      mem_bytes(TYPE_MAT,0,sizeof(MAT));
      mem_numvar(TYPE_MAT,1);
   }
   
   matrix->m = m;       matrix->n = matrix->max_n = n;
   matrix->max_m = m;   matrix->max_size = m*n;
#ifndef SEGMENTED
   if ((matrix->base = NEW_A(m*n,Real)) == (Real *)NULL )
   {
      free(matrix);
      error(E_MEM,"m_get");
   }
   else if (mem_info_is_on()) {
      mem_bytes(TYPE_MAT,0,m*n*sizeof(Real));
   }
#else
   matrix->base = (Real *)NULL;
#endif
   if ((matrix->me = (Real **)calloc(m,sizeof(Real *))) == 
       (Real **)NULL )
   {    free(matrix->base); free(matrix);
    error(E_MEM,"m_get");
     }
   else if (mem_info_is_on()) {
      mem_bytes(TYPE_MAT,0,m*sizeof(Real *));
   }
   
#ifndef SEGMENTED
   /* set up pointers */
   for ( i=0; i<m; i++ )
     matrix->me[i] = &(matrix->base[i*n]);
#else
   for ( i = 0; i < m; i++ )
     if ( (matrix->me[i]=NEW_A(n,Real)) == (Real *)NULL )
       error(E_MEM,"m_get");
     else if (mem_info_is_on()) {
    mem_bytes(TYPE_MAT,0,n*sizeof(Real));
       }
#endif
   
   return (matrix);
}

Why do they first allocate base and use it to allocate me ? Also, if you've read the source code, please tell me the use of SEGMENTED in this library. The declaration is in configure file.

The matrix structure is defined in matrix.h and m_get() is in memory.c.

#define NEW_A(num,type) ((type *)calloc((size_t)(num),sizeof(type))) — kneelb4darth, Feb 15 '16 at 17:07

score 2 · Accepted Answer · edited May 23 '17 at 11:45

If I read the code correctly, me is a one-dimensional array of pointers to the starts of the "rows" in the matrix: matrix->me[i] = &(matrix->base[i*n]); fills the array with the addresses of the respective row at index i.

That makes it possible to use me with two indices as if it were indeed a two-dimensional array, e.g. write double d = myMatrix.me[row][column];. Of course sizeof(myMatrix.me) or taking the address of me's elements behaves differently than for a proper 2-dimensional array.

If SEGMENTED is defined, the memory in the matrix will be allocated row by row and will not be contiguous (and base will be null). If you want to pass a raw 2-dimensional array to a library unaware of the MAT type the elements will need to lie contiguous in memory. If you just index normally it doesn't matter.

As an aside (I know that this is code from a library, but still): It's funny that just today I was directed to an interesting solution to this problem by Jens Gustedt, https://stackoverflow.com/a/12805980/3150802. He uses a C99 feature, variable length arrays. Interestingly he doesn't actually create any variable length array (instead, he just mallocs like your matrix class) but just uses the run-time length in the declaration. Very interesting. The A in his post plays the role of me in your matrix class, i.e. it's a pointer to one-dimensional arrays which can be indexed twice to obtain a single number from the matrix.

The trick is that it *looks* like one, but is completely 1-dimensional. Scrutinize this expression `&(matrix->base[i*n])`: It is a single address value, a number if you want, and it is assigned to an element of `me`. (Using a *single* index). It's true: At the memory location defined by this address there is a whole (one-dimensional) array of REAL again, namely a row of the matrix. Indexing `me` in `me[row]` yields that row, so that we can apply an index to the result of the first index in order to get one element of that row, as in `(me[row])[column]`. — Peter - Reinstate Monica, Feb 15 '16 at 17:29

difficulty understanding the data structure in meschach library

1 Answers1