How to efficiently store a triangular matrix in memory?

Question

I want to store a lower triangular matrix in memory, without storing all the zeros. The way I have implemented it is by allocating space for i + 1 elements on the ith row. However, I am new to dynamic memory allocation in C and something seems to be wrong with my first allocation.

int main ()
{
    int i, j;
    int **mat1;
    int dim;

    scanf("%d", &dim);
    *mat1 = (int**) calloc(dim, sizeof(int*));

    for(i = 0; i < dim; i++)
    mat1[i] = (int*) calloc(i + 1, sizeof(int));

    for(i = 0; i < dim; i++)
    {
        for(j = 0; j < i + 1; j++)
        {
            scanf("%d", &mat1[i][j]);
        }
    }

    /* Print the matrix without the zeros*/
    for(i = 0; i < dim; i++)
    {
        for(j = 0; j < (i + 1); j++)
        {
            printf("%d%c", mat1[i][j], j != (dim-1) ? ' ' : '\n');
        }
    }

    return 0;
}

Answer 1

If you want to conserve space and the overhead of allocating every row of the matrix, you could implement a triangular matrix by using clever indexing of a single array.

A lower triangular matrix (including diagonals) has the following properties:

Dimension   Matrix    Elements/row   Total elements
1           x . . .   1              1
2           x x . .   2              3
3           x x x .   3              6
4           x x x x   4              10
...

The total number of elements for a given dimension is:

size(d) = 1 + 2 + 3 + ... + d  =  (d+1)(d/2)

If you lay the rows out consecutively in a single array, you can use the formula above to calculate the offset of a given row and column (both zero-based) inside the matrix:

index(r,c) = size(r-1) + c

The formulas above are for the lower triangular matrix. You can access the upper matrix as if it was a lower matrix by simply reversing the indexes:

index((d-1)-r, (d-1)-c)

If you have concerns about changing the orientation of the array, you can devise a different offset calculation for the upper array, such as:

uindex(r,c) = size(d)-size(d-r) + c-r

Sample code:

#include <time.h>
#include <stdio.h>
#include <stdlib.h>

#define TRM_SIZE(dim) (((dim)*(dim+1))/2)
#define TRM_OFFSET(r,c) (TRM_SIZE((r)-1)+(c))
#define TRM_INDEX(m,r,c) ((r)<(c) ? 0 : (m)[TRM_OFFSET((r),(c))])
#define TRM_UINDEX(m,r,c,d) ((r)>(c)?0:(m)[TRM_SIZE(d)-TRM_SIZE((d)-(r))+(c)-(r)])
#define UMACRO 0


int main (void)
{
  int i, j, k, dimension;
  int *ml, *mu, *mr;

  printf ("Enter dimension: ");
  if (!scanf ("%2d", &dimension)) {
    return 1;
  }

  ml = calloc (TRM_SIZE(dimension), sizeof *ml);
  mu = calloc (TRM_SIZE(dimension), sizeof *mu);
  mr = calloc (dimension*dimension, sizeof *mr);
  if (!ml || !mu || !mr) {
    free (ml);
    free (mu);
    free (mr);
    return 2;
  }

  /* Initialization */

  srand (time (0));
  for (i = 0; i < TRM_SIZE(dimension); i++) {
    ml[i] = 100.0*rand() / RAND_MAX;
    mu[i] = 100.0*rand() / RAND_MAX;
  }

  /* Multiplication */

  for (i = 0; i < dimension; i++) {
    for (j = 0; j < dimension; j++) {
      for (k = 0; k < dimension; k++) {
        mr[i*dimension + j] +=
#if UMACRO
          TRM_INDEX(ml, i, k) *
          TRM_UINDEX(mu, k, j, dimension);
#else
          TRM_INDEX(ml, i, k) *
          TRM_INDEX(mu, dimension-1-k, dimension-1-j);
#endif
      }
    }
  }

  /* Output */

  puts ("Lower array");
  for (i = 0; i < dimension; i++) {
    for (j = 0; j < dimension; j++) {
      printf (" %2d", TRM_INDEX(ml, i, j));
    }
    putchar ('\n');
  }
  puts ("Upper array");
  for (i = 0; i < dimension; i++) {
    for (j = 0; j < dimension; j++) {
#if UMACRO
      printf (" %2d", TRM_UINDEX(mu, i, j, dimension));
#else
      printf (" %2d", TRM_INDEX(mu, dimension-1-i, dimension-1-j));
#endif
    }
    putchar ('\n');
  }
  puts ("Result");
  for (i = 0; i < dimension; i++) {
    for (j = 0; j < dimension; j++) {
      printf (" %5d", mr[i*dimension + j]);
    }
    putchar ('\n');
  }

  free (mu);
  free (ml);
  free (mr);

  return 0;
}

Note that this is a trivial example. You could extend it to wrap the matrix pointer inside a structure that also stores the type of the matrix (upper or lower triangular, or square) and the dimensions, and write access functions that operate appropriately depending on the type of matrix.

For any non-trivial use of matrices, you should probably use a third-party library that specializes in matrices.

Answer 2

mat1 = calloc(dim,sizeof(int*));

mat1 is a double pointer.You need to allocate memory for your array of pointers and later you need to allocate memory to each of your pointers individually.No need to cast calloc()

Answer 3

You are dereferencing mat1 at line 8 before it has even been set to point anywhere. You are allocating an array of pointers to int, but you are not assigning that to mat1 but to the dereference of mat1, which is uninitialized, we don't know what it points to.

So this line:

// ERROR: You are saying an unknown memory location should have the value of calloc.
*mat1 = (int**)calloc(dim,sizeof(int*));

Should change to:

// OK: Now you are assigning the allocation to the pointer variable.
mat1 = (int**)calloc(dim,sizeof(int*));