BLAS库如何直接与升压多阵列一起使用?

问题描述

给出两个矩阵AB,并给出C的相乘结果。

#include "boost/multi_array.hpp"
typedef boost::multi_array<double,2> matrix;
int m=5;
int n=6;
int k=7;

matrix A(boost::extents[m][k]);
matrix B(boost::extents[k][n]);
matrix C(boost::extents[m][n]);

如何调用dgemm库中的blas函数来计算AB的矩阵乘积? 我知道uBLAS部分boost库,armadilloMTL 4eigen和其他一些为blas函数提供方便包装的库。这里的问题是如何直接在多阵列上调用dgemm

解决方法

您可以访问连续元素存储。

原型是

void cblas_dgemm( CBLAS_LAYOUT layout,CBLAS_TRANSPOSE TransA,CBLAS_TRANSPOSE TransB,const int M,const int N,const int K,const double alpha,const double *A,const int lda,const double *B,const int ldb,const double beta,double *C,const int ldc )

因此,让我们填写它:

cblas_dgemm(
    CBLAS_LAYOUT::CblasRowMajor,CBLAS_TRANSPOSE::CblasNoTrans,m,n,k,1.0,// alpha
    A.data(),A.shape()[1],B.data(),B.shape()[1],0.0,// beta
    C.data(),C.shape()[1]);

那是使用docs here来获得有关alpha / beta的一些指导:

演示

#include <boost/multi_array.hpp>
typedef boost::multi_array<double,2> matrix;

#include <iostream>
namespace io { // for debug output
    auto& dump(std::ostream& os,double v) { return os << v; }

    template <typename R> auto& dump(std::ostream& os,R const& r) {
        std::string_view sep = "";
        os << "{";
        for (auto const& el : r) { dump(os << sep,el); sep = ","; }
        return os << "}";
    }
}

std::ostream& operator<<(std::ostream& os,matrix const& m) { return io::dump(os,m); }

// demo
#include <cblas-netlib.h>
#include <numeric> // iota

int main() {
    constexpr auto m=5,n=6,k=7;

    matrix A(boost::extents[m][k]);
    matrix B(boost::extents[k][n]);
    matrix C(boost::extents[m][n]);

    std::iota(A.data(),A.data() + A.num_elements(),0);
    std::iota(B.data(),B.data() + B.num_elements(),50);
    std::iota(C.data(),C.data() + C.num_elements(),100);

    std::cout << "A: " << A << "\nB: " << B << "\n";
    assert(A.storage_order().all_dims_ascending());

    /*
     * void cblas_dgemm(
     *   CBLAS_LAYOUT layout,*   CBLAS_TRANSPOSE TransA,*   CBLAS_TRANSPOSE TransB,*   const int M,*   const double alpha,*   const double *A,*   const double *B,*   const double beta,*   double *C,const int ldc )
     */

    cblas_dgemm(
        CBLAS_LAYOUT::CblasRowMajor,// alpha
        A.data(),// beta
        C.data(),C.shape()[1]);

    std::cout << "C:\n" << C << "\n";
}

哪些印刷品:

A: {{0,1,2,3,4,5,6},{7,8,9,10,11,12,13},{14,15,16,17,18,19,20},{21,22,23,24,25,26,27},{28,29,30,31,32,33,34}}
B: {{50,51,52,53,54,55},{56,57,58,59,60,61},{62,63,64,65,66,67},{68,69,70,71,72,73},{74,75,76,77,78,79},{80,81,82,83,84,85},{86,87,88,89,90,91}}
C:
{{1596,1617,1638,1659,1680,1701},{4928,4998,5068,5138,5208,5278},{8260,8379,8498,8617,8736,8855},{11592,11760,11928,12096,12264,12432},{14924,15141,15358,15575,15792,16009}}

此操作以Wolfram Alpha签出:

enter image description here