Matrix multiplication algorithm problem (C ++)

Hi, I am helping my friend write a program that multiplies two upper triangular matrices without matrix expansion. When I say no matrix expansion, I mean not filling the bottom of the top triangular matrix with zeros (the goal is to preserve space). I am reading in a matrix from a file that only has upper triangular values ​​for matrices with the lower part omitted.

What I'm trying to accomplish is to write a function that has given an array and a pair of indices, will return the element that the expanded matrix will have at that location (0 for diagonals below and values ​​above diagonal). I think about it using the usual matrix multiplication algorithm that uses this function whenever it needs to access an element. I've been working on this for several hours and I can't think of a way to convert the indices of a double matrix (i, j) (i goes along the rows) to indices of a single array and vice versa (as a reminder I'm using a one-dimensional array to store an upper triangular matrix) ... Any help would be greatly appreciated!

+3


source to share


3 answers


// mat is an array containing the upper triangle data for a square matrix of size n
// returns element at (i,j), or 0 for the lower triangle
int getFromTriangle(const int* mat, int n, int i, int j)
{
  if (i > j) {
    return 0; // lower triangle
  } else {
    return mat[j + (i*n) - i*(i+1)/2];
  }
}

      

The proposal if

takes care of the bottom triangle. The clause else

calculates the array index as follows:



j + (i*n) - i*(i+1)/2

      

This is just the usual function of the index of a rectangular matrix minus the offset, which is exactly the i

th triangular number , since in any row i

, the store has omitted triangle(i)

elements.

0


source


You can split the algorithm into smaller, easily understandable chunks.



// Get the array index given the rank, the row, and the column.
int getArrayIndex(int rank, int row, int col)
{
   return (row*rank + col - col*(col+1)/2);
}

// Get the term of a matrix, given the rank, the row, and the column.
int getMatrixTerm(int a[], int rank, int row, int col)
{
   if ( col < row )
   {
      return 0;
   }
   else
   {
      return a[getArrayIndex(rank, row, col)];
   }
}

// Get the term for a row and column resulting from mulitiplication.
int getMultipliedTerm(int a[], int b[], int rank, int row, int col)
{
   int term = 0;
   int k = 0;
   for ( ; k < rank; ++k )
   {
      term += getMatrixTerm(a, rank, row, k)*getMatrixTerm(b, rank, k, col);
   }
   return term;
}

// Set the term in c given the rank, the row, and the column.
void setMultipliedTerm(int a[], int b[], int c[], int rank, int row, int col)
{
   if ( j >= i )
   {
      c[getArrayIndex(rank, i, j)] = getMultipliedTerm(a, b, rank, i, j);
   }
}

// High level function to multiply two upper triangular matrices
// The lower part of the matrix is not stored at all.
void multiply(int a[], int b[], int c[], int rank)
{
   int i = 0;
   int j = 0;
   for ( i = 0; i < rank; ++i )
   {
      for ( j = 0; j < rank; ++j )
      {
         setMultipliedTerm(a, b, c, rank, i, j);
      }
   }
}

      

0


source


if p = & a [0] [0] then a [i] [j] is the same as * (p + i * row_size + j) in C ++. where p is a pointer to the same data type as the matrix elements. Hope this is what you want and that you know pointers. all the best.

0


source







All Articles