矩阵访问和CPU矩阵乘法优化

4

我正在使用JNI帮助下开发一些优化的Java矩阵封装。需要确认的是,您能给出一些有关矩阵优化的提示吗?我要实现的是:

矩阵可以表示为四组缓冲区/数组,一个用于水平访问,一个用于垂直访问,一个用于对角线访问,以及一个命令缓冲区,仅在需要时计算矩阵元素。这里是一个示例。

Matrix signature: 

0  1  2  3  

4  5  6  7

8  9  1  3

3  5  2  9

First(hroizontal) set: 
horSet[0]={0,1,2,3} horSet[1]={4,5,6,7} horSet[2]={8,9,1,3} horSet[3]={3,5,2,9}

Second(vertical) set:
verSet[0]={0,4,8,3} verSet[1]={1,5,9,5} verSet[2]={2,6,1,2} verSet[3]={3,7,3,9}

Third(optional) a diagonal set:
diagS={0,5,1,9} //just in case some calculation needs this

Fourth(calcuation list, in a "one calculation one data" fashion) set:
calc={0,2,1,3,2,5} --->0 means multiply by the next element
                       1 means add the next element
                       2 means divide by the next element
                       so this list means
                       ( (a[i]*2)+3 ) / 5  when only a[i] is needed.
Example for fourth set: 
A.mult(2),   A.sum(3),  A.div(5), A.mult(B)
(to list)   (to list)  (to list) (calculate *+/ just in time when A is needed )
 so only one memory access for four operations.
 loop start
 a[i] = b[i] * ( ( a[i]*2) +3 ) / 5  only for A.mult(B)
 loop end

如上所述,当需要访问列元素时,第二个集合提供连续访问,没有跳跃。水平访问可以通过第一个集合实现相同的效果。

这样做会使一些事情变得更容易,一些事情变得更困难:

 Easier: 
 **Matrix transpozing operation. 
 Just swapping the pointers horSet[x] and verSet[x] is enough.

 **Matrix * Matrix multiplication.
 One matrix gives one of its horizontal set and other matrix gives vertical buffer.
 Dot product of these must be highly parallelizable for intrinsics/multithreading.
 If the multiplication order is inverse, then horizontal and verticals are switched.

 **Matrix * vector multiplication.
 Same as above, just a vector can be taken as horizontal or vertical freely.

 Harder:
 ** Doubling memory requirement is bad for many cases.
 ** Initializing a matrix takes longer.
 ** When a matrix is multiplied from left, needs an update vertical-->horizontal
 sets if its going to be multiplied from right after.(same for opposite)
 (if a tranposition is taken between, this does not count)


 Neutral:
 ** Same matrix can be multiplied with two other matrices to get two different
 results such as A=A*B(saved in horizontal sets)   A=C*A(saved in vertical sets)
 then A=A*A gives   A*B*C*A(in horizontal) and C*A*A*B (in vertical) without
 copying A. 

 ** If a matrix always multiplied from left or always from right, every access
 and multiplication will not need update and be contiguous on ram.

 ** Only using horizontals before transpozing, only using verticals after, 
 should not break any rules.

主要目的是拥有一个大小为(8的倍数,8的倍数)的矩阵,并使用多个线程应用AVX指令集(每个线程同时处理一组)。
我只实现了向量*向量点积。如果您编程方面很精通,我会进一步学习。
我编写的点积(使用Intrinsics)比展开循环版本快6倍(后者是一一相乘的两倍),但在包装器中启用多线程时卡在存储带宽上限处(8x --> 使用近20GB/s,接近DDR3的极限)。已经尝试过OpenCL,对于CPU来说速度较慢,但对于GPU来说非常好。
谢谢。
编辑: "块矩阵"缓冲区如何执行?在乘法大矩阵时,小补丁以特殊方式相乘,并且可能使用缓存以减少主内存访问。但这需要在竖直-水平-对角和该块之间的矩阵乘法之间进行更多更新。

1
获得高度优化的矩阵代数代码的一种方法是使用表达式模板,在不执行每个特定步骤的情况下将整个内容压缩到专门的方法中。 - Pixelchemist
使用子矩阵和一些线性代数来利用缓存更多(块算法),怎么样?这对于垂直和水平缓冲区方法可行吗? - huseyin tugrul buyukisik
2个回答

1

一些库使用表达式模板来实现针对级联矩阵操作的非常具体且优化的函数应用。

C++编程语言也有一个关于“融合操作”的简短章节(第4版29.5.4)。

这使得可以像这样连接语句:

M = A*B.transp(); // where M, A, B are matrices

在这种情况下,您需要有3个类:

class Matrix;

class Transposed
{
public:
  Transposed(Matrix &matrix) : m_matrix(matrix) {}
  Matrix & obj (void) { return m_matrix; }
private:
  Matrix & m_matrix;
};

class MatrixMatrixMulTransPosed
{
public:
  MatrixMatrixMulTransPosed(Matrix &matrix, Transposed &trans) 
    : m_matrix(matrix), m_transposed(trans.obj()) {}
  Matrix & matrix (void) { return m_matrix; }
  Matrix & transposed (void) { return m_transposed; }
private:
  Matrix & m_matrix;
  Matrix & m_transposed;
};

class Matrix
{
  public:
    MatrixMatrixMulTransPosed operator* (Transposed &rhs)
    { 
      return MatrixMatrixMulTransPosed(*this, rhs); 
    }

    Matrix& operator= (MatrixMatrixMulTransPosed &mmtrans)
    {
      // Actual computation goes here and is stored in this.
      // using mmtrans.matrix() and mmtrans.transposed()
    }
};

你可以将这个概念推进,使其能够针对任何关键计算拥有专门的功能。

我正在测试向量x向量,它在45毫秒内(1.4GHz fx8150)将两个6400万元素的向量相乘。您是否也期望矩阵乘法具有类似的性能?(例如将大小为256x512和512 x 512的两个矩阵相乘) - huseyin tugrul buyukisik
我不是线性代数大师,但我猜“普通”的矩阵乘法在跨步数据访问的情况下会遭受很多缓存未命中,但你显然愿意通过第二个数据集来解决这个问题。然而,你应该小心,避免引入任何人为的瓶颈,因为两种表示需要匹配,并且你需要有更新例程。 - Pixelchemist
好的,第一次乘法已经完成。 在1.4GHz fx8150的单个核心上,使用1024x1024矩阵,乘法花费了1.3秒。 你认为这是否有任何缓存未命中? - huseyin tugrul buyukisik
当我将问题规模增加到4096x4096(工作量增加了64倍),将CPU频率增加到4.0GHz,并启用C#的Parallel.For并行化时(是的,我将代码迁移到了具有易于并行化功能的C#),乘法运算花费了14.9秒。 - huseyin tugrul buyukisik

1
这实际上相当于缓存转置。听起来你打算急切地这样做;我建议只在需要时计算转置,并记住它以防再次需要。这样,如果你从未需要它,则不会进行计算。

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接