矩阵与向量高性能计算实现与代码范例

为了方便用户进行线性代数与矢量分析，金字塔提供了一整套封装的高度矢量优化线性代数方面函数，用户可以不借助Python也可以在PEL中进行简易矩阵及向量计算。金字塔会根据用户CPU自动启用SSE,AVX,AVX-512指令集，以及高度并行化执行，可以发挥CPU最大利用率，执行效率大大提升。

使用事项说明
      1、该类函数的向量及矩阵的输入参数只能是VARIABLE申明的数组类型，使用前应注意转换。
      2、函数计算需依赖科学计算组件，该组件您需要单独安装部署，使用前请先进行在线安装。
      安装步骤 ： PEL公式编辑器界面->插入菜单->安装扩展科学计算组件。
      3、该类库函数还在不断扩充中，如果你有其他的功能函数需求，欢迎及时在本贴中提出，以便我们进行扩充安排。

1、Mfsum(X,N)
计算向量x的元素大小之和
函数实现功能： y=∑x
参数维度说明： x为向量(n*1)(矩阵也可以看做向量)

1. VARIABLE:x[]=0;
   
2.      x:=c;
   
3.      rets:Mfsum(x,5);

复制代码

2、Mfxpy(x,y,alpha,n)
计算向量与标量的积，然后加到结果上
函数实现功能：y:=a*x+y
x : 向量
y : 向量(输出参数）
    alpha : 标量
    n : 指定向量x,y的元素个数
    返回值 : 最终更新到向量y

1. alpha:=1;
   
2. VARIABLE : x[5]=(0,1,2,3,4);
   
3. VARIABLE : y[5]=(1,2,3,4,5);
   
4. If IsLastbar Then
   
5. Begin
   
6.         MFXPY(x,y,ALPHA,5);
   
7.         for i = 1 to 5 do
   
8.         mSGout(1,numtostr(y[i],0));//结果将输出 1,3,5,7,9
   
9.     End   

复制代码

3、Mfdot(x,y,n)
计算向量-向量的点乘
x : 向量
y : 向量
n：向量x,y数量

1. VARIABLE : x[5]=(0,1,2,3,4);
   
2. VARIABLE : y[5]=(1,2,3,4,5);
   
3. RETS:MFDOT(x,y,5);//结果将输出40

复制代码

4、Mfnrm2(x,n)
计算一个向量的欧几里得范数(Euclidean norm)
函数实现功能：res = ||x||
     x:向量
     n:向量x个数

1. VARIABLE : x[5]=(1,2,3,4,5);
   
2. RETS:MFNRM2(x,5); //结果将输出 7.416

复制代码

5、Mfscal(x,a,n)
计算向量和标量的乘积
函数实现功能：x=a*x
x : 向量（输入/输出参数）
a : 标量a
n: 向量x个数
返回值：最终更新到向量x

1. VARIABLE : x[5]=(1,2,3,4,5);
   
2. a:=0.1;
   
3. If IsLastbar Then
   
4. Begin
   
5.         MFSCAL(x,a,5);
   
6.         for i = 1 to 5 do
   
7.       Msgout(1,numtostr(x[i],1));  //结果将输出 0.1,0.2,0.3,0.4,0.5
   
8. End  

复制代码

6、Mfswap(x,y,n)
交换向量值, 完成向量x和y的交换
x : 向量
y : 向量
n：向量x,y数量

1. VARIABLE : x[5]=(0,1,2,3,4);
   
2. VARIABLE : y[5]=(1,2,3,4,5);
   
3. If IsLastbar Then
   
4. Begin
   
5.     Mfswap(x,y,5);
   
6.     for i = 1 to 5 do
   
7.     msgout(1,numtostr(y[i],0)); //输出交换后的y数组
   
8. End  

复制代码

7、mfamax(x,n)
找到绝对值最大的元素的索引
x: 向量
n: 向量x个数
返回值：基于1索引的向量元素位置

1. VARIABLE : x[5]=(0,1,2,3,4);
   
2. Rets:MFAMAX(x,5); //结果将返回5

复制代码

8、Mfamin(x,n)
找到绝对值最小的元素的索引
x: 向量
n: 向量x个数
返回值：基于1索引的向量元素位置
例：见Mfamax函数
9、mfAdd(a,b,r,n)
向量元素相加
a : 向量
b : 向量
r：向量(输出值)
n: 向量a,b,r个数

1. VARIABLE : a[5]=(1,2,3,4,5);
   
2. VARIABLE : b[5]=(1,2,3,4,5);
   
3. VARIABLE : r[5]=0;
   
4. If IsLastbar Then
   
5. Begin
   
6.     MFADD(a,b,r,5);
   
7.     for i = 1 to 5 do
   
8.         msgout(1,numtostr(r[i],0)); //结果将返回2，4，6，8，10
   
9. End  

复制代码

10、MfSub(a,b,r,n)
向量元素相减
a : 向量;b : 向量;r : 向量(输出值);n: 向量a,b,r个数
例：见mfAdd函数

11、MfMul(a,b,r,n)
向量元素相乘
a : 向量;b : 向量;r：向量(输出值);n: 向量a,b,r个数
例：见mfAdd函数

12、Mfgemv(a,x,y,t,m,n,alpha,beta)
矩阵与向量的乘积
函数实现功能：
t =0 y := alpha*A*x + beta*y,
t =1 y := alpha*A'*x + beta*y,
t =2 y := alpha*conjg(A')*x + beta*y,
a :  m*n的矩阵
x : 向量;y : 向量(输入和输出参数)
m : 矩阵A的行数;n : 矩阵A的列数（向量x维度）
alpha : 标量;beta : 标量
返回值：最终更新到向量y

1. m:=2;//矩阵行数
   
2. n:=5;//矩阵列数
   
3. alpha:=1;
   
4. mbeta:=1;
   
5. VARIABLE : a[m*n]=( 0,1,2,3,4,5,6,7,8,9 ); //定义一个2行5列矩阵
   
6. VARIABLE : x[n]=(1,2,3,4,5); //向量数组
   
7. VARIABLE : y[m]=0;
   
8. If IsLastbar Then
   
9.       Begin
   
10.       //2*5的矩阵与5*1的向量相乘
    
11.           MFGEMV(a,x,y,0,m,n,alpha,mbeta);
    
12.           for i = 1 to m do
    
13.           msgout(1,numtostr(y[i],0)); //结果将输出40, 115
    
14. End

复制代码

Mfger(x,y,a,m,n,alpha)
矩阵的一阶更新（矩阵-向量乘积）
函数实现功能：A := alpha*x*y'+ A,
x : 向量(m行)
y : 向量(n列)
a : m×n 的矩阵(输入和输出参数)
m : 矩阵A的行数
n : 矩阵A的列数
alpha : 标量
返回值：最终更新到矩阵a

1. m:=2;//矩阵行数
   
2. n:=5;//矩阵列数
   
3. alpha:=1;
   
4. VARIABLE : x[m]=(1,2);
   
5. VARIABLE : y[n]=(2,3,4,5,6);
   
6. VARIABLE : a[m*n]=( 1,1,1,1,1,1,1,1,1,1 ); //定义一个2行5列矩阵
   
7. If IsLastbar Then
   
8.        Begin
   
9.        //5*1向量乘以1*5向量，加上5*5矩阵
   
10.        MFGER(x,y,a,m,n,alpha);
    
11.        for i = 1 to m*n do
    
12.        msgout(1,numtostr(a[i],0));
    
13. End  
    
14. //结果将输出矩阵为
    
15. 3  4  5  6  7
    
16.      5  7  9 11 13

复制代码

Mfgemm(A,B,C,alpha,beta,m,n,k)
矩阵与矩阵乘积
函数实现功能：C := alpha*op(A)*op(B) + beta*C,
A: 矩阵
B: 矩阵
C:矩阵（输入和输出参数）
alpha: 标量
beta : 标量
m：矩阵A和矩阵C的行数
n：矩阵B和矩阵C的列数
k：矩阵A的列数和矩阵B的行数
返回值：最终更新到矩阵C

1. m:=2;//矩阵A和矩阵C的行数
   
2. n:=3;//矩阵B和矩阵C的列数
   
3. k:=2;//矩阵A的列数和矩阵B的行数
   
4. alpha:=1;
   
5. betas:=1;
   
6. VARIABLE : a[m*n]=(1,2,3,  4,5,6); //矩阵维度 2*3
   
7. VARIABLE : b[n*k]=(1,1,  1,1,  1,1); //矩阵维度 3*2(计算时候转置)
   
8. VARIABLE : x[m*k]=(0,1,  2,3);      //矩阵维度 2*2
   
9. If IsLastbar Then
   
10.        Begin
    
11.        MFGEMM(a,b,x,alpha,betas,m,n,k);
    
12.        for i = 1 to m*k do
    
13.        msgout(1,numtostr(x[i],0));
    
14. End  
    
15. //结果将输出矩阵为
    
16. 6 7
    
17. 1718

复制代码

Mffft(X,Y,N,T)
快速傅立叶变换
X: 输入源向量数组
Y: 输出向量数组，以复数方式返回(由2个数字组成)，因此数组申明大小必须为X的2倍
N: X向量变换源数组大小
T:0向前变换 1向后变换

1. n:=16; //16点信号变换
   
2. VARIABLE : x[n]=(15, 32, 9, 222, 118, 151, 5, 7, 56, 233, 56, 121, 235, 89, 98, 111);
   
3. VARIABLE : y[n*2]=0;  //以复数方式返回结果，因此需要2倍容量
   
4. If IsLastbar Then
   
5. Begin
   
6.         MFFFT(x,y,n,0);
   
7.         for i = 1 to n*2 do
   
8.         msgout(1,numtostr(y[i],3));
   
9. End  

复制代码

矩阵和多维数组的使用介绍
矩阵也是多维数组的一种(二维数组)，矩阵有行和列两项，在内存中的排列是连续的，所以一维数组也可以进行2维数组矩阵进行表达。
访问矩阵的计算公式为 ： (访问行-1)*总列数+访问列

范例：

1. m:=2; //矩阵行数
   
2. n:=5; //矩阵列数
   
3. VARIABLE:x[m*n]=(1,2,3,4,5,   6,7,8,9,10);//定义并初始化2行5列矩阵
   
4. //取2行3列数据，二维数组表示[2][3];
   
5. y:=x[(2-1)*n+3]; //y=8
   
6. s:y;
   
7. //取1行5列数据，二维数组表示[1][5]
   
8. y:=x[(1-1)*n+5]; //y=5
   
9. t:y;

复制代码