Intel MKL基础(4)MKL函数举例(BLAS and Sparse BLAS)

时间:2022-05-03 16:04:41
MKL实现了BLAS和Spare BLAS和BLAS-like的函数。这部分函数主要分为以下部分:

•BLAS Level 1 Routines (vector-vector operations) 
•BLAS Level 2 Routines (matrix-vector operations) 
•BLAS Level 3 Routines (matrix-matrix operations) 
•Sparse BLAS Level 1 Routines (vector-vector operations). 
•Sparse BLAS Level 2 and Level 3 Routines (matrix-vector and matrix-matrix operations) 
•BLAS-like Extensions 

一、BLAS

(1)BLAS函数命名规则

<character> <name> <mod> ( )

其中:

character:

表示数据类型,一般为:s实数,单精度;c复数,双精度;d实数,单精度;z复数,双精度。有些函数会对character组合,比如sc、dz等等,组合的时候,第一个表示返回值,第二个字母表示参数。

name:

对于BLAS level 1: 表示操作类型,比如?dot,?rot,?swap等函数。(说明:MKL手册中?开头的函数,这个问号就取代character,表示这一系列的函数。)

对于BLAS level 2和3: name反应了矩阵的参数类型,主要有ge(general matrix)、gb(general band matrix带状矩阵)、sy(symmetric matrix对称式矩阵)、sp(symmetric matrix (packed storage))、sb(symmetric band matrix)、he(Hermitian matrix埃尔米特矩阵,自共轭矩阵)、hp(Hermitian matrix (packed storage))、hb(Hermitian band matrix)、tr(triangular matrix三角矩阵)、tp(triangular matrix (packed storage))、tb(triangular band matrix)等。

mod:

可选,提供操作的额外信息。不同的level的函数可取的值不一样,参考手册查看更多信息。

总结:可见MKL函数的命名规则是很清晰的。

(2)矩阵的保存方式

矩阵,在MKL中可以用一维或者二维的数组保存。
Full storage:用二维数组保存矩阵,矩阵的元素aij保存在数据元素的a(i,j)中。
Packed storage:用一维数组保存。
band storage for band matrices:用二维矩阵保存,带状矩阵。

(3)编程接口

MKL手册中发现函数说明都是Fortran的规范,比如BLAS里面的函数,事实上,MKL提供了和Fortran一样的函数接口。

C语言中,使用BLAS需要包含的头文件是:#include <mkl_blas.h>,其中的函数和Fortran接口是一致的。

为了使得接口更符合C语言的格式,MKL提供了一个C语言格式的函数接口,即CBLAS,需要包含#include <mkl_cblas.h>,其提供的接口函数名为Fortran函数名接口加上前缀"cblas_"的形式,但是实际上其参数返回值等可能会有一些不同,是为了更符合C语言的编程习惯。具体参考下一部分的例子即可知道,应该尽量采用C语言的接口形式,更容易理解和使用。

(3)BLAS Level 1举例

BLAS Level 1的函数用于进行向量-向量的操作。比如求向量元素的绝对值之和、向量相加、向量交换、向量元素最大值和最小值等等操作。

下面以?asum函数为例来了解:

res = sasum(n, x, incx)
res = scasum(n, x, incx)
res = dasum(n, x, incx)
res = dzasum(n, x, incx)

其中,不同的函数,其?部分不一样,表示其传递的参数类型不一样,上面已经分析过命名规则。这个函数用于对向量中的一组元素的绝对值求和,如果是复数参数,则其元素会计算实部和虚部的绝对值并求和。即res = |Re x(1)| + |Im x(1)| + |Re  x(2)| + |Im  x(2)|+ ... + |Re  x(n)| + |Im x(n)|。其中,参数n表示要参与计算的元素的个数,x为向量起始地址,incx为对向量的递增取值间隔,注意,使用时要自行保证不要出现越界。另外,上面的函数中n和incx都要求是指针类型,所以不符合C语言的编程习惯,可以使用C语言接口,从而符合习惯。参考下面的例子:

#include "stdafx.h"
#include <stdio.h>
#include <malloc.h>
#include <mkl_blas.h>
#include <mkl_cblas.h>

#define N 100
void initVector(float* v)
{
for(int i = 0;i < N;i++)
v[i] = -i*1.0f;
}

int main(int argc, char* argv[])
{
float* vector = (float*) malloc(sizeof(float)*N);
initVector(vector);

int a1 = N;
int incx1 = 1;
float ret1 = sasum(&a1, vector, &incx1); // FORTRAN 77接口形式
printf("Result of sasum: %lf\n", ret1);

int a2 = N/2;
int incx2 = 2;
float ret2 = sasum(&a2, vector, &incx2); // FORTRAN 77接口形式
printf("Result of sasum: %lf\n", ret2);

float ret3 = cblas_sasum(N, vector, 1); //C接口形式
printf("Result of sasum: %lf\n", ret3);

free(vector);

return 0;
}
这个例子,是计算一个单精度浮点数向量的。

下面是计算复数的例子:

#include "stdafx.h"
#include <stdio.h>
#include <malloc.h>
#include <mkl_cblas.h>

#define N 10
void initVector(MKL_Complex16* v)
{
for(int i = 0;i < N;i++)
{
v[i].real = -i*1.0f;
v[i].imag = i*1.0f;
}
}

int main(int argc, char* argv[])
{
MKL_Complex16* vector = (MKL_Complex16*) malloc(sizeof(MKL_Complex16)*N);
initVector(vector);

double ret3 = cblas_dzasum(N, vector, 1); // 复数
printf("Result of sasum: %lf\n", ret3);

free(vector);

return 0;
}

还有其它一些BLAS Level 1的函数,都是用于处理vector-vecotr相关的BLAS问题的,更多内容参考MKL手册。

(4)BLAS Level 2

BLAS Level 2的函数,主要用于矩阵-向量(matrix-vector)操作。

(5)BLAS Level 3

BLAS Level 2的函数,主要用于矩阵-矩阵(matrix-matrix)操作。

二、Sparse BLAS

(1)命名规则

Level1 和Level2/3有不同的命名规则,详细参考手册。

(2)Sparse BLAS Level 1

主要用于针对压缩形式保存的稀疏向量的相关向量操作。

(3) Sparse BLAS Level 2 and Level 3

Level 2 主要用于进行稀疏矩阵和密集向量之间的操作。

Level 3 主要用于稀疏矩阵和密集矩阵之间的操作。

三、BLAS-like Extensions

作为扩展,提供一些扩展的函数,进行一些数据操作等相关的函数,比如简单的向量加法得到一个新的向量等。


总结:这里主要是简单的了解了一下MKL手册中的BLAS and Sparse BLAS的基本内容,主要是举例BLAS Level1的例子,了解MKL函数的特点。另外,这些函数,不仅仅在MKL中有提供,这些都是“标准”的一些线性代数相关的函数,本身BLAS、LAPACK、ScaLAPACK等等的这些名称就是一些库的名称,都是在高性能领域很有名的库,MKL主要是对其进行了优化,所以,这些函数的接口都是和“标准”一致的。