矩阵乘法算法问题（C++）

Question

你好，我正在帮助我的朋友编写一个程序，用于不扩展矩阵的情况下两个上三角矩阵相乘。在说不扩展矩阵时，我是指不用零填充上三角矩阵的下部（目标是节省空间）。我从文件中读取矩阵，该文件仅包含省略了矩阵下部分的上三角值。

我试图编写一个函数，该函数给定一个数组和一对索引，将返回扩展矩阵在该位置的元素（下三角为0，上三角为值）。我考虑使用通常的矩阵乘法算法，每当需要访问元素时使用此函数。我已经花费了几个小时的时间，但我不能想出一种将双倍矩阵索引（i，j）（i沿行）转换为单个数组索引和反之亦然的方法（提醒一下，我使用单维数组存储上三角矩阵）。非常感谢你的帮助！

- JustAnotherUser

---

对于 j == i，主对角线上的元素位于 0、n、n + (n-1)、n + (n-1) + (n-2) 等位置（计算该序列的第 i 个元素留给读者作为练习）。对于 j > i，只需将 j - i 添加到前面计算的结果中即可。 - Igor Tandetnik

Answer 1

// 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];
  }
}

if语句负责处理下三角部分。而else语句则按照以下方式计算数组索引：

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

这只是常规的矩形索引函数减去一个偏移量，该偏移量恰好是第个三角形数，因为在任何一行中，存储都省略了triangle(i)个元素。

- John Zwinck

Answer 2

你可以将算法分解成更小、更易理解的部分。

// 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);
      }
   }
}

- R Sahu

Answer 3

如果p = &a [0] [0]，则a [i] [j]在c ++中等同于*(p + i * row_size + j)。 其中p是指向与矩阵元素相同类型的数据的指针。 希望这就是你想要的，并且你了解指针。 一切顺利。

- pii_ke