Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How to efficiently store a triangular matrix in memory?

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;
}
like image 568
rares985 Avatar asked Dec 09 '14 13:12

rares985


People also ask

How is an upper triangular matrix stored in the computer memory?

When an upper-triangular matrix is stored in upper-triangular-packed storage mode, the upper triangle of the matrix is stored, including the diagonal, in a one-dimensional array. The upper triangle is packed by columns.

How do you store a lower triangular matrix?

Such a matrix can be implemented in an optimized manner. The efficient way to store the lower triangular matrix of size N: Count of non-zero elements = 1 + 2 + 3 + … + N = N * (N + 1) /2.

What is memory representation of a triangular matrix?

The four storage modes used for storing triangular matrices are described in the following: Upper-Triangular-Packed Storage Mode. Lower-Triangular-Packed Storage Mode. Upper-Triangular Storage Mode. Lower-Triangular Storage Mode.


2 Answers

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.

like image 63
Nisse Engström Avatar answered Oct 18 '22 10:10

Nisse Engström


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()

like image 45
Gopi Avatar answered Oct 18 '22 10:10

Gopi