为什么增强矩阵乘法比我的慢？

Question

我已经实现了一个矩阵乘法boost::numeric::ublas::matrix(参见我的完整工作增强代码)

Result result = read ();

boost::numeric::ublas::matrix<int> C;
C = boost::numeric::ublas::prod(result.A, result.B);

另一个使用标准算法(参见完整的标准代码):

vector< vector<int> > ijkalgorithm(vector< vector<int> > A, 
                                    vector< vector<int> > B) {
    int n = A.size();

    // initialise C with 0s
    vector<int> tmp(n, 0);
    vector< vector<int> > C(n, tmp);

    for (int i = 0; i < n; i++) {
        for (int k = 0; k < n; k++) {
            for (int j = 0; j < n; j++) {
                C[i][j] += A[i][k] * B[k][j];
            }
        }
    }
    return C;
}

这就是我测试速度的方法:

time boostImplementation.out > boostResult.txt
diff boostResult.txt correctResult.txt

time simpleImplementation.out > simpleResult.txt
diff simpleResult.txt correctResult.txt

两个程序都读取了一个包含两个2000 x 2000矩阵的硬编码文本文件.这两个程序都是用这些标志编译的:

g++ -std=c++98 -Wall -O3 -g $(PROBLEM).cpp -o $(PROBLEM).out -pedantic

我15秒对我实施和在4分钟用于升压的实现!

编辑:编译后

g++ -std=c++98 -Wall -pedantic -O3 -D NDEBUG -DBOOST_UBLAS_NDEBUG library-boost.cpp -o library-boost.out

ikj算法得到28.19 秒,Boost得到60.99秒.所以Boost仍然相当慢.

为什么提升比我的实现慢得多？

Answer 1

如TJD所指出的那样,可以通过调试后者的功能来部分地解释uBLAS版本的较慢性能.

以下是uBLAS版本调试时间:

real    0m19.966s
user    0m19.809s
sys     0m0.112s

以下是关闭调试的uBLAS版本所花费的时间(-DNDEBUG -DBOOST_UBLAS_NDEBUG添加了编译器标志):

real    0m7.061s
user    0m6.936s
sys     0m0.096s

因此,关闭调试,uBLAS版本几乎快3倍.

剩余的性能差异可以通过引用uBLAS FAQ的以下部分来解释"为什么uBLAS比(atlas-)BLAS慢得多":

ublas的一个重要设计目标是尽可能通用.

这种普遍性几乎总是带来成本.特别地,prod函数模板可以处理不同类型的矩阵,例如稀疏矩阵或三角矩阵.幸运的是,uBLAS提供了针对密集矩阵乘法优化的替代方案,特别是axpy_prod和block_prod.以下是比较不同方法的结果:

ijkalgorithm   prod   axpy_prod  block_prod
   1.335       7.061    1.330       1.278

正如您所看到的那样,axpy_prod并且block_prod比您的实现更快.只测量没有I/O的计算时间,去除不必要的复制并仔细选择块大小block_prod(我使用64)可以使差异更加深刻.

另请参阅uBLAS常见问题解答和有效的uBlas以及常规代码优化.

Answer 2

我相信,你的编译器没有足够的优化.uBLAS代码大量使用模板和模板需要大量使用优化.我在发布模式下通过MS VC 7.1编译器为1000x1000矩阵运行代码,它给了我

10.064用于uBLAS

7.851s为矢量

差异仍然存在,但绝不是压倒性的.uBLAS的核心概念是延迟评估,因此prod(A, B)仅在需要时评估结果,例如prod(A, B)(10,100)将立即执行,因为实际上只会计算一个元素.因此,实际上没有可以优化的整个矩阵乘法的专用算法(见下文).但你可以帮助图书馆一点点,宣布

matrix<int, column_major> B;

将4.426单手操作系统的运行时间减少到与你的功能相媲美的时间.这个声明使得在乘法矩阵时更加顺序地访问内存,从而优化缓存使用.

PS已经阅读了uBLAS文档到最后;),您应该已经发现实际上有一个专用函数可以立即乘以整个矩阵.2个功能 - axpy_prod和opb_prod.所以

opb_prod(A, B, C, true);

即使在未经优化的row_major B矩阵中,也可以在几秒内执行,8.091并且与矢量算法相同

PPS还有更多优化:

C = block_prod<matrix<int>, 1024>(A, B);

4.4无论B是column_还是row_ major,都在s中执行.考虑描述:"函数block_prod是为大密集矩阵设计的." 为特定任务选择特定工具!