openblas 矩阵算法库简介

1 向量与向量

2 矩阵与向量

3 矩阵与矩阵

openblas 是一个开源的矩阵计算库，包含了诸多的精度和形式的矩阵计算算法。就精度而言，包括float和double，两种数据类型的数据，其矩阵调用函数也是不一样。不同矩阵，其计算方式也是有所不同，（姑且认为向量也是一维矩阵），例如，向量与向量之间的计算，向量与矩阵之间的计算，矩阵与矩阵之间的计算。

按照不同的计算需求openblas有如下几种

1 向量与向量

对于不同的精度，使用的函数也是不一样。例如cblasdaxpy表示double类型数据，cblasdaxpy表示float类型数据

头文件

"cblas.h"

描述

cblas_?axpy系列函数执行vector-vector的操作，即：

y := alpha*x + y;

其中alpha是标量，x，y是向量。

参数说明

n Specifies the number of elements in vectors x and y.

a Specifies the scalar a.

x Array, size at least (1 + (n-1)*abs(incx)).

incx Specifies the increment for the elements of x. 一般设置为1

y Array, size at least (1 + (n-1)*abs(incy)).

incy Specifies the increment for the elements of y. 一般设置为1

2 矩阵与向量

矩阵与向量的乘法种类也有很多，常规也就是上述三种，以cblas_sgemv为例：

头文件：

"cblas.h"

描述：

cblas_?ger系列函数执行matrix-vector的操作，即：

y := alpha*A*x + beta*y,

其中alpha，beta是标量，x，y是向量（x(m*1)，y(n*1)），A是一个m*n的矩阵。

详细接口：

void cblas_sgemv (const CBLAS_LAYOUT Layout, const CBLAS_TRANSPOSE trans, const MKL_INT m, const MKL_INT n, const float alpha, const float *a, const MKL_INT lda, const float *x, const MKL_INT incx, const float beta, float *y, const MKL_INT incy)

详细接口的具体含义：

Layout：布局，也就是指定矩阵的排列方式，按照行优先(CblasRowMajor)还是列优先(CblasColMajor)，需要指出的是，在函数使用中，虽然是向量与矩阵的乘法，但是实际上，传递的是两个指针向量。这样对于多维矩阵而言，指明其排列方式，十分有必要。参数（CblasRowMajor，CblasColMajor）

trans：转置特定操作：

if trans=CblasNoTrans, then y := alphaAx + beta * y;

if trans=CblasTrans, then y := alphaA'x + beta * y;

if trans=CblasConjTrans, then y := alpha conjg(A')x + beta*y.

m：表示矩阵A的行m

n：表示矩阵A的列n

a：表示一个矩阵，size = lda *k

For Layout = CblasColMajor, k is n，

For Layout = CblasRowMajor, k is m.

alpha：一个标量，通常为1

x：矩阵向量

incx，incy，一般为1

beta: 一个标量参数，通常设置为 0，之前试过设置为1，直接导致了结果不正确；在矩阵或是向量元素的值较大的时候，beta这个值，为0还是1影响不大；

示例程序1：

void matrix_vector(){float array[6] = { 1,2,3,4,5,6 };float x[3] = { 1,2,3 };int m = 2; int n = 3;float y[2] = {};cblas_sgemv(CblasRowMajor, CblasNoTrans, m, n, alpha, array, n, x, 1, beta, y, inc_y);//cblas_sgemv(CblasColMajor, CblasTrans, n, m, alpha, array, n, x, 1, beta, y, inc_y);//两个结果一致
}//打印结果
14,32

示例程序2：

void matrix_vector()
{int i = 0;double x[2] = { 1.0,2.0 };double y[3] = { 2.0,1.0,3.0};double A[6] = { 0 };blasint rows = 2, cols = 3;double alpha = 10;blasint inc_x = 1, inc_y = 1;blasint lda = 2;//矩阵按列优先存储//A <== alpha*x*y' + A （y'表示y的转置）cblas_dger(CblasColMajor, rows, cols, alpha, x, inc_x, y, inc_y, A, lda);for (i = 0; i<6; i++)std::cout << A[i] << " ";std::cout << "\n";
}
//打印结果
20 40 10 20 30 60

3 矩阵与矩阵

矩阵与矩阵的乘积，相较于前两种较为复杂，但是用途却比较大，上面能算的，它行；上面不能算的，它也行；

头文件：

"cblas.h"

描述：

C := alphaop(A)op(B) + betaC,

其中：op（A）表示A矩阵或是转置等
矩阵
B是一个KN矩阵
C是一个MN 矩阵

A是一个M

接口描述：

void cblas_sgemm(OPENBLAS_CONST enum CBLAS_ORDER Order, OPENBLAS_CONST enum CBLAS_TRANSPOSE TransA, OPENBLAS_CONST enum CBLAS_TRANSPOSE TransB, OPENBLAS_CONST blasint M, OPENBLAS_CONST blasint N, OPENBLAS_CONST blasint K,

OPENBLAS_CONST float alpha, OPENBLAS_CONST float *A, OPENBLAS_CONST blasint lda, OPENBLAS_CONST float *B, OPENBLAS_CONST blasint ldb, OPENBLAS_CONST float beta, float *C, OPENBLAS_CONST blasint ldc);

详细参数描述：

cblas_sgemm(order,transA,transB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC);

详细参数与上述的矩阵与向量乘法类似

以此函数为例：

cblas_sgemm(order,transA,transB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC);

alpha =1,beta =0 的情况下，等于两个矩阵相成。

第一参数 oreder 候选值有ClasRowMajow 和ClasColMajow 这两个参数决定一维数组怎样存储在内存中,

一般用ClasRowMajow

参数 transA和transB ：表示矩阵A，B是否进行转置。候选参数 CblasTrans 和CblasNoTrans.

参数M：表示 A或C的行数。如果A转置，则表示转置后的行数

参数N：表示 B或C的列数。如果B转置，则表示转置后的列数。

参数K：表示 A的列数或B的行数（A的列数=B的行数）。如果A转置，则表示转置后的列数。

参数LDA：表示A的列数，与转置与否无关。

参数LDB：表示B的列数，与转置与否无关。

参数LDC：始终=N

NOTE: 首先判定是行优先还是列优先；

再次,依照lda，ldb分别将两个矩阵按行、列分开

再次,判定是否进行转置操作

输出结果；

示例程序1：

#include <vector>
#include <iostream>
#include "caffe/util/math_functions.hpp"
using namespace std;
int main() {const int M=4;const int N=2;const int K=3;const float alpha=1;const float beta=0;const int lda=M;const int ldb=K;const int ldc=N;const float A[K*M]={1,2,3,4,5,6,7,8,9,8,7,6};const float B[N*K]={5,4,3,2,1,0};float C[M*N];cblas_sgemm(CblasRowMajor, CblasTrans, CblasTrans, M, N, K, alpha, A, lda, B, ldb, beta, C, ldc);for(int i=0;i<M;i++){for(int j=0;j<N;j++){cout<<C[i*N+j]<<" ";}cout<<endl;}
」
输出：
52 7
58 10
64 13
70 16

示例程序2：

#include <vector>
#include <iostream>
#include "caffe/util/math_functions.hpp"
using namespace std;
int main() {const int M=4;const int N=2;const int K=3;const float alpha=1;const float beta=0;const int lda=K;const int ldb=N;const int ldc=N;const float A[M*K]={1,2,3,4,5,6,7,8,9,8,7,6};const float B[K*N]={5,4,3,2,1,0};float C[M*N];cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, M, N, K, alpha, A, lda, B, ldb, beta, C, ldc);for(int i=0;i<M;i++){for(int j=0;j<N;j++){cout<<C[i*N+j]<<" ";}   cout<<endl;}
} 结果：
14 8
41 26
68 44
67 46