第26部分- Linux ARM汇编浮点和向量

最新推荐文章于 2021-06-04 11:31:07 发布

badman250

最新推荐文章于 2021-06-04 11:31:07 发布

阅读量3.7k

点赞数

分类专栏： 64位 ARM处理器汇编技术系列

本文链接：https://blog.csdn.net/notbaron/article/details/106578635

版权

64位 ARM处理器汇编技术系列专栏收录该内容

29 篇文章 14 订阅

订阅专栏

第26部分- Linux ARM汇编浮点和向量

二进制浮点数是一个实数的近似表示，由三个部分组成：符号，尾数和指数。

1.01110 是20 + 2-2 + 2-3 + 2-4 = 1.43750(10)

VFPv2支持两个IEEE 754数字：Binary32和Binary64，通常以其C类型分别为float和double或单精度和双精度。

在单精度浮点中，尾数为23位（归一化数字为整数1的+1），指数为8位（因此指数范围为-126到127）。

在双精度浮点中，尾数为52位（+1），指数为11位（因此，指数范围为-1022至1023）。

32位的时候浮点计算是以协处理器的方式存在，64位的时候并不是协处理器存在了。

VFP和NEON在ARMV8-A中是一个整体。

C代码

#include <stdio.h>
#include <string.h>
#include <math.h>
#define M 4
#define N 4
int main(int argc, char const *argv[])
{
	float a[M][N] = {{0.1,0.2,0.0,0.1},{0.2,0.1,0.3,0.0},{0.0,0.3,0.1,0.5},{0.0,0.6,0.4,0.1}};
	float b[N][M] = {{4.92,2.54,-0.63,-1.75},{3.02,-1.51,-0.87,1.35},{-4.29,2.14,0.71,0.71},{-0.95,0.48,2.38,-0.95}};
	int re = 0;
	float c[M][M] = {{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0}};
	for (int i = 0; i < M; ++i)
	{
		for (int j = 0; j < M; ++j)
		{
			for (int k = 0; k < N; ++k)
			{
				c[i][j] += a[i][k] * b[k] [j];
			}
		}
	}
	for (int i = 0; i < M; ++i)
	{
		for (int j = 0; j < M; ++j)
		{
			printf("%5.3f ", c[i][j]);
			if (j == M - 1)
			{
				printf("\n");
			}
		}
	}
	return 0;
}

# gcc -o matmulc matmul.c -std=c99

$./matmulc

1.001 -0.000 0.001 0.000

-0.001 0.999 0.000 -0.002

0.002 0.001 1.000 0.001

0.001 -0.002 0.000 0.999

32位矩阵乘

v = <v0, v1, … vr-1>

w = <w0, w1, …, wr-1>

矩阵：

v · w = v0×w0 + v1×w1 + … + vr-1×wr-1

C代码中矩阵相乘片段如下：

float A[N][M]; // N rows of M columns each row
float B[M][P]; // M rows of P columns each row
// Result
float C[N][P];
 
for (int i = 0; i < N; i++) // for each row of the result
{
  for (int j = 0; j < P; j++) // and for each column
  {
    C[i][j] = 0; // Initialize to zero
    // Now make the dot matrix of the row by the column
    for (int k = 0; k < M; k++)
       C[i][j] += A[i][k] * B[k][j];
  }
}

看下汇编的矩阵乘

原始代码

3个循环，i,j,k。

/* -- matmul.s */
.data
mat_A: .float 0.1, 0.2, 0.0, 0.1
       .float 0.2, 0.1, 0.3, 0.0
       .float 0.0, 0.3, 0.1, 0.5 
       .float 0.0, 0.6, 0.4, 0.1
mat_B: .float  4.92,  2.54, -0.63, -1.75
       .float  3.02, -1.51, -0.87,  1.35
       .float -4.29,  2.14,  0.71,  0.71
       .float -0.95,  0.48,  2.38, -0.95
mat_C: .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
 
format_result : .asciz "Matrix result is:\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n"
 
.text
naive_matmul_4x4:
    /* r0 address of A
       r1 address of B
       r2 address of C
    */
    push {r4, r5, r6, r7, r8, lr} /* 保存整型寄存器*/
    /* First zero 16 single floating point */
    /* In IEEE 754, all bits cleared means 0 */
    mov r4, r2;//矩阵C地址给r4
    mov r5, #16;//矩阵C条目数量
    mov r6, #0;//矩阵C都赋值为0
    b .Lloop_init_test//初始化矩阵C
    .Lloop_init :
      str r6, [r4], +#4   /* *r4 ← r6 then r4 ← r4 + 4 */
    .Lloop_init_test:
      subs r5, r5, #1
      bge .Lloop_init
 
    /* 使用r4为i，r5为j，r6 为k*/
    mov r4, #0 /* r4 ← 0 ,i=0*/
    .Lloop_i:  /* 循环i */
      cmp r4, #4  /* if r4 == 4 循环i结束 */
      beq .Lend_loop_i
      mov r5, #0  /* r5 ← 0,j从0开始 */
      .Lloop_j: /*循环j*/
       cmp r5, #4 /* if r5 == 4循环j结束 */
        beq .Lend_loop_j
        /* 计算C[i][j]地址并加载到s0 */
        /* C[i][j] 地址是 C + 4*(4 * i + j) */
        mov r7, r5               /* r7 ← r5. r7=j */
        adds r7, r7, r4, LSL #2  /* r7 ← r7 + (r4 << 2). 即r7 ← j + i * 4. We multiply i by the row size (4 elements) */
        adds r7, r2, r7, LSL #2  /* r7 ← r2 + (r7 << 2).即r7 ← C + 4*(j + i * 4)， 用原始大小乘以(j + i * 4) */
        vldr s0, [r7] /* s0 ← *r7 ，元素加载到s0*/
        mov r6, #0 /* r6 ← 0 */
        .Lloop_k :  /* 循环 k */
          cmp r6, #4 /* if r6 == 4循环 k 结束*/
          beq .Lend_loop_k
          /* 计算a[i][k]地址，并加载到s1*/
          /* a[i][k] 地址a + 4*(4 * i + k) */
          mov r8, r6  /* r8 ← r6. r8=k */
          adds r8, r8, r4, LSL #2  /* r8 ← r8 + (r4 << 2). 即r8 ← k + i * 4 */
          adds r8, r0, r8, LSL #2  /* r8 ← r0 + (r8 << 2). 即r8 ← a + 4*(k + i * 4) */
          vldr s1, [r8]    /* s1 ← *r8 ，即a[i][j]放到s1中*/
          /* 计算b[k][j]地址，并加载到s2 */
          /* b[k][j] 地址b + 4*(4 * k + j) */
          mov r8, r5               /* r8 ← r5. 即r8 ← j */
          adds r8, r8, r6, LSL #2  /* r8 ← r8 + (r6 << 2). 即r8 ← j + k * 4 */
          adds r8, r1, r8, LSL #2  /* r8 ← r1 + (r8 << 2). 即r8 ← b + 4*(j + k * 4) */
          vldr s2, [r8]  /* s1 ← *r8 ，b[k][j]放到s2*/
          vmul.f32 s3, s1, s2      /* s3 ← s1 * s2 */
          vadd.f32 s0, s0, s3      /* s0 ← s0 + s3 */
 
          add r6, r6, #1           /* r6 ← r6 + 1,k++ */
          b .Lloop_k               /* 循环k下个迭代 */
        .Lend_loop_k: /* Here ends loop k */
        vstr s0, [r7]      /* 保存s0 回 C[i][j] */
        add r5, r5, #1  /* r5 ← r5 + 1,j++ */
        b .Lloop_j /*下个j迭代*/
       .Lend_loop_j: /* Here ends loop j */
       add r4, r4, #1 /* r4 ← r4 + 1，i++*/
       b .Lloop_i     /*下个i迭代 */
    .Lend_loop_i: /* Here ends loop i */
    pop {r4, r5, r6, r7, r8, lr}  /* Restore integer registers */
    bx lr /* Leave function */
 
.globl main
main:
    push {r4, r5, r6, lr}  /*保存整型寄存器*/
    vpush {d0-d1}          /* 保存浮点寄存器*/
 
    /* 准备调用函数naive_matmul_4x4 */
    ldr r0, addr_mat_A  /* r0 ← a */
    ldr r1, addr_mat_B  /* r1 ← b */
    ldr r2, addr_mat_C  /* r2 ← c */
    bl naive_matmul_4x4
 
    /* 打印结果矩阵*/
    ldr r4, addr_mat_C  /* r4 ← c矩阵C地址 */
    vldr s0, [r4] /* s0 ← *r4. 第一个参数s0 ← c[0][0] */
    vcvt.f64.f32 d1, s0 /* 转成双精度d1 ← s0*/
    vmov r2, r3, d1      /* {r2,r3} ← d1 */
 
    mov r6, sp     /* 保存sp,r6 ← sp */
    mov r5, #1  /* 从1-15，共60个字节 */
    add r4, r4, #60 /* 矩阵c最后一个地址是c[3][3] */
    .Lloop:
        vldr s0, [r4] /* s0 ← *r4. 加载当前项*/
        vcvt.f64.f32 d1, s0 /*转成双精度d1 ← s0*/
        sub sp, sp, #8      /* 在堆栈预留一个双精度的空间*/
        vstr d1, [sp]       /* 将当前项保存到堆栈中 */
        sub r4, r4, #4      /* 矩阵地址减4个字节到前一项*/
        add r5, r5, #1      /* 待处理的矩阵项数目*/
        cmp r5, #16         /* 判断是会否结束*/
        bne .Lloop
 
    ldr r0, addr_format_result /* 格式化字符串*/
    bl printf /* 调用函数printf */
    mov sp, r6  /* 恢复sp寄存器*/
 
    mov r0, #0
    vpop {d0-d1}
    pop {r4, r5, r6, lr}
    bx lr
 
addr_mat_A : .word mat_A
addr_mat_B : .word mat_B
addr_mat_C : .word mat_C
addr_format_result : .word format_result

as -g -o matmul.o matmul.s -mfpu=vfpv2

gcc -o matmul matmul.o

输出如下：

$./matmul

Matrix result is:

1.00 -0.00 0.00 0.00

-0.00 1.00 0.00 -0.00

0.00 0.00 1.00 0.00

0.00 -0.00 0.00 1.00

向量优化

由于A [i] [k]和A [i] [k + 1]是内存中的连续元素，因此访问A [i] [k]有多重加载。可以完全避免所有循环k并执行从A [i] [0]到A [i] [3]的4元素加载。

通过这种方法，可以一次删除4个操作，从而完全删除循环k。必须修改fpscr，以便将len字段设置为4（并在离开函数时将其恢复为1）。

ARMv7 使用NEON instructions (也叫做高级SIMD)

编译的时候加入-mfpu=neon

.data
mat_A: .float 0.1, 0.2, 0.0, 0.1
       .float 0.2, 0.1, 0.3, 0.0
       .float 0.0, 0.3, 0.1, 0.5 
       .float 0.0, 0.6, 0.4, 0.1
mat_B: .float  4.92,  2.54, -0.63, -1.75
       .float  3.02, -1.51, -0.87,  1.35
       .float -4.29,  2.14,  0.71,  0.71
       .float -0.95,  0.48,  2.38, -0.95
mat_C: .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
 
format_result : .asciz "Matrix result is:\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n"
 
.text
naive_vec_matmul_4x4:
    /* r0 address of A
       r1 address of B
       r2 address of C
    */
push {r4, r5, r6, r7, r8, lr} /* 保存整型寄存器*/
    vpush {s16-s19}  //保存s16-s19
    vpush {s24-s27}

    mov r4, r2;//矩阵C地址给r4
    mov r5, #16;//矩阵C条目数量
    mov r6, #0;//矩阵C都赋值为0
    b .Lloop_init_test//初始化矩阵C
    .Lloop_init :
      str r6, [r4], +#4   /* *r4 ← r6 then r4 ← r4 + 4 */
    .Lloop_init_test:
      subs r5, r5, #1
      bge .Lloop_init
 
    /* 使用r4为i，r5为j，r6 为k*/
    mov r4, #0 /* r4 ← 0 ,i=0*/
    .Lloop_i:  /* 循环i */
      cmp r4, #4  /* if r4 == 4 循环i结束 */
      beq .Lend_loop_i
      mov r5, #0  /* r5 ← 0,j从0开始 */
      .Lloop_j: /*循环j*/
       cmp r5, #4 /* if r5 == 4循环j结束 */
        beq .Lend_loop_j
        /* 计算C[i][j]地址并加载到s0 */
        /* C[i][j] 地址是 C + 4*(4 * i + j) */
        mov r7, r5               /* r7 ← r5. r7=j */
        adds r7, r7, r4, LSL #2  /* r7 ← r7 + (r4 << 2). 即r7 ← j + i * 4. We multiply i by the row size (4 elements) */
        adds r7, r2, r7, LSL #2  /* r7 ← r2 + (r7 << 2).即r7 ← C + 4*(j + i * 4)， 用原始大小乘以(j + i * 4) */

        /* a[i][0]地址 */
        mov r8, r4, LSL #2
        adds r8, r0, r8, LSL #2
        //vldmia r8, {s8-s11}  /* 矢量加载{s8,s9,s10,s11} ← {a[i][0], a[i][1], a[i][2], a[i][3]} */
        vldmia r8, {s8,s9,s10,s11}  /* 矢量加载{s8,s9,s10,s11} ← {a[i][0], a[i][1], a[i][2], a[i][3]} */
 
        /* b[0][j]地址 */
        mov r8, r5               /* r8 ← r5. This is r8 ← j */
        adds r8, r1, r8, LSL #2  /* r8 ← r1 + (r8 << 2). This is r8 ← b + 4*(j) */
        vldr s16, [r8]       /* s16 ← *r8. 即s16 ← b[0][j] */
        vldr s17, [r8, #16] /* s17 ← *(r8 + 16). 即s17 ← b[1][j] */
        vldr s18, [r8, #32] /* s18 ← *(r8 + 32). 即s17 ← b[2][j] */
        vldr s19, [r8, #48] /* s19 ← *(r8 + 48). 即s17 ← b[3][j] */

	vmul.f32 s0,s8,s16
	vmla.f32 s0,s9,s17
	vmla.f32 s0,s10,s18
	vmla.f32 s0,s11,s19
 
        vstr s0, [r7]            /* 保存s0 到 C[i][j] */

        add r5, r5, #1  /* r5 ← r5 + 1,j++ */
        b .Lloop_j /*下个j迭代*/
       .Lend_loop_j: /* Here ends loop j */
       add r4, r4, #1 /* r4 ← r4 + 1，i++*/
       b .Lloop_i     /*下个i迭代 */
.Lend_loop_i: /* Here ends loop i */

    vpop {s24-s27}                
    vpop {s16-s19}

    pop {r4, r5, r6, r7, r8, lr}  /* Restore integer registers */
    bx lr /* Leave function */
 
.globl main
main:
    push {r4, r5, r6, lr}  /*保存整型寄存器*/
    vpush {d0-d1}          /* 保存浮点寄存器*/
 
    /* 准备调用函数naive_matmul_4x4 */
    ldr r0, addr_mat_A  /* r0 ← a */
    ldr r1, addr_mat_B  /* r1 ← b */
    ldr r2, addr_mat_C  /* r2 ← c */
    bl naive_vec_matmul_4x4
 
    /* 打印结果矩阵*/
    ldr r4, addr_mat_C  /* r4 ← c矩阵C地址 */
    vldr s0, [r4] /* s0 ← *r4. 第一个参数s0 ← c[0][0] */
    vcvt.f64.f32 d1, s0 /* 转成双精度d1 ← s0*/
    vmov r2, r3, d1      /* {r2,r3} ← d1 */
 
    mov r6, sp     /* 保存sp,r6 ← sp */
    mov r5, #1  /* 从1-15，共60个字节 */
    add r4, r4, #60 /* 矩阵c最后一个地址是c[3][3] */
    .Lloop:
        vldr s0, [r4] /* s0 ← *r4. 加载当前项*/
        vcvt.f64.f32 d1, s0 /*转成双精度d1 ← s0*/
        sub sp, sp, #8      /* 在堆栈预留一个双精度的空间*/
        vstr d1, [sp]       /* 将当前项保存到堆栈中 */
        sub r4, r4, #4      /* 矩阵地址减4个字节到前一项*/
        add r5, r5, #1      /* 待处理的矩阵项数目*/
        cmp r5, #16         /* 判断是会否结束*/
        bne .Lloop
 
    ldr r0, addr_format_result /* 格式化字符串*/
    bl printf /* 调用函数printf */
    mov sp, r6  /* 恢复sp寄存器*/
 
    mov r0, #0
    vpop {d0-d1}
    pop {r4, r5, r6, lr}
    bx lr
 
addr_mat_A : .word mat_A
addr_mat_B : .word mat_B
addr_mat_C : .word mat_C
addr_format_result : .word format_result

as -g -o matmul_vec.o matmul_vec.s -mfpu=neon

gcc -o matmul_vec matmul_vec.o

这里使用 neon向量，vmul，还没优化。

寄存器并发

如果同时处理C [i] [j]和C [i] [j + 1]可以填充每个库中的所有8个寄存器。

.data
mat_A: .float 0.1, 0.2, 0.0, 0.1
       .float 0.2, 0.1, 0.3, 0.0
       .float 0.0, 0.3, 0.1, 0.5 
       .float 0.0, 0.6, 0.4, 0.1
mat_B: .float  4.92,  2.54, -0.63, -1.75
       .float  3.02, -1.51, -0.87,  1.35
       .float -4.29,  2.14,  0.71,  0.71
       .float -0.95,  0.48,  2.38, -0.95
mat_C: .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
 
format_result : .asciz "Matrix result is:\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n"
 
.text
naive_vec_matmul_4x4:
    /* r0 address of A
       r1 address of B
       r2 address of C
    */
push {r4, r5, r6, r7, r8, lr} /* 保存整型寄存器*/
    vpush {s16-s19}  //保存s16-s19
    vpush {s24-s27}

    mov r4, r2;//矩阵C地址给r4
    mov r5, #16;//矩阵C条目数量
    mov r6, #0;//矩阵C都赋值为0
    b .Lloop_init_test//初始化矩阵C
    .Lloop_init :
      str r6, [r4], +#4   /* *r4 ← r6 then r4 ← r4 + 4 */
    .Lloop_init_test:
      subs r5, r5, #1
      bge .Lloop_init
 
    /* 使用r4为i，r5为j，r6 为k*/
    mov r4, #0 /* r4 ← 0 ,i=0*/
    .Lloop_i:  /* 循环i */
      cmp r4, #4  /* if r4 == 4 循环i结束 */
      beq .Lend_loop_i
      mov r5, #0  /* r5 ← 0,j从0开始 */
      .Lloop_j: /*循环j*/
       cmp r5, #4 /* if r5 == 4循环j结束 */
        beq .Lend_loop_j
        /* 计算C[i][j]地址并加载到s0 */
        /* C[i][j] 地址是 C + 4*(4 * i + j) */
        mov r7, r5               /* r7 ← r5. r7=j */
        adds r7, r7, r4, LSL #2  /* r7 ← r7 + (r4 << 2). 即r7 ← j + i * 4. We multiply i by the row size (4 elements) */
        adds r7, r2, r7, LSL #2  /* r7 ← r2 + (r7 << 2).即r7 ← C + 4*(j + i * 4)， 用原始大小乘以(j + i * 4) */

        /* a[i][0]地址 */
        mov r8, r4, LSL #2
        adds r8, r0, r8, LSL #2
        vldmia r8, {s8,s9,s10,s11}  /* 矢量加载{s8,s9,s10,s11} ← {a[i][0], a[i][1], a[i][2], a[i][3]} */
 
        /* b[0][j]地址 */
        mov r8, r5               /* r8 ← r5. This is r8 ← j */
        adds r8, r1, r8, LSL #2  /* r8 ← r1 + (r8 << 2). This is r8 ← b + 4*(j) */
        vldr s16, [r8]       /* s16 ← *r8. 即s16 ← b[0][j] */
        vldr s17, [r8, #16] /* s17 ← *(r8 + 16). 即s17 ← b[1][j] */
        vldr s18, [r8, #32] /* s18 ← *(r8 + 32). 即s17 ← b[2][j] */
        vldr s19, [r8, #48] /* s19 ← *(r8 + 48). 即s17 ← b[3][j] */

       /* b[0][j-1]地址 */
        add r8, r5,#1               /* r8 ← r5+1. This is r8 ← j+1 */
        adds r8, r1, r8, LSL #2  /* r8 ← r1 + (r8 << 2). This is r8 ← b + 4*(j) */
        vldr s20, [r8]       /* s20 ← *r8. 即s20 ← b[0][j] */
        vldr s21, [r8, #16] /* s21 ← *(r8 + 16). 即s21 ← b[1][j] */
        vldr s22, [r8, #32] /* s22 ← *(r8 + 32). 即s22 ← b[2][j] */
        vldr s23, [r8, #48] /* s23 ← *(r8 + 48). 即s23 ← b[3][j] */

	vmul.f32 s0,s8,s16
	vmla.f32 s0,s9,s17
	vmla.f32 s0,s10,s18
	vmla.f32 s0,s11,s19

        vmul.f32 s1,s8,s20
	vmla.f32 s1,s9,s21
	vmla.f32 s1,s10,s22
	vmla.f32 s1,s11,s23

        vstmia r7,{s0-s1}      /* 保存s0 到 C[i][j],C[i][j] */

        add r5, r5, #2  /* r5 ← r5 + 1,j++ */
        b .Lloop_j /*下个j迭代*/
       .Lend_loop_j: /* Here ends loop j */
       add r4, r4, #1 /* r4 ← r4 + 1，i++*/
       b .Lloop_i     /*下个i迭代 */
.Lend_loop_i: /* Here ends loop i */

    vpop {s24-s27}                
    vpop {s16-s19}

    pop {r4, r5, r6, r7, r8, lr}  /* Restore integer registers */
    bx lr /* Leave function */
 
.globl main
main:
    push {r4, r5, r6, lr}  /*保存整型寄存器*/
    vpush {d0-d1}          /* 保存浮点寄存器*/
 
    /* 准备调用函数naive_matmul_4x4 */
    ldr r0, addr_mat_A  /* r0 ← a */
    ldr r1, addr_mat_B  /* r1 ← b */
    ldr r2, addr_mat_C  /* r2 ← c */
    bl naive_vec_matmul_4x4
 
    /* 打印结果矩阵*/
    ldr r4, addr_mat_C  /* r4 ← c矩阵C地址 */
    vldr s0, [r4] /* s0 ← *r4. 第一个参数s0 ← c[0][0] */
    vcvt.f64.f32 d1, s0 /* 转成双精度d1 ← s0*/
    vmov r2, r3, d1      /* {r2,r3} ← d1 */
 
    mov r6, sp     /* 保存sp,r6 ← sp */
    mov r5, #1  /* 从1-15，共60个字节 */
    add r4, r4, #60 /* 矩阵c最后一个地址是c[3][3] */
    .Lloop:
        vldr s0, [r4] /* s0 ← *r4. 加载当前项*/
        vcvt.f64.f32 d1, s0 /*转成双精度d1 ← s0*/
        sub sp, sp, #8      /* 在堆栈预留一个双精度的空间*/
        vstr d1, [sp]       /* 将当前项保存到堆栈中 */
        sub r4, r4, #4      /* 矩阵地址减4个字节到前一项*/
        add r5, r5, #1      /* 待处理的矩阵项数目*/
        cmp r5, #16         /* 判断是会否结束*/
        bne .Lloop
 
    ldr r0, addr_format_result /* 格式化字符串*/
    bl printf /* 调用函数printf */
    mov sp, r6  /* 恢复sp寄存器*/
 
    mov r0, #0
    vpop {d0-d1}
    pop {r4, r5, r6, lr}
    bx lr
 
addr_mat_A : .word mat_A
addr_mat_B : .word mat_B
addr_mat_C : .word mat_C
addr_format_result : .word format_result

as -g -o matmul_vec2.o matmul_vec2.s -mfpu=neon

gcc -o matmul_vec2 matmul_vec2.o

访问优化

修改外层循环的次序，将K层移动到最外层，可以到的如下：

float A[N][N];
float B[M][N];
// Result
float C[N][N];
 
for (int i = 0; i < N; i++)
  for (int j = 0; j < N; j++)
    C[i][j] = 0;
 
for (int k = 0; k < N; k++)
  for (int i = 0; i < N; i++)
    for (int j = 0; j < N; j++)
       C[i][j] += A[i][k] * B[k][j];

最里面的循环就是

for (int k = 0; k < N; k++)
  for (int i = 0; i < N; i++)
  {
     C[i][0] += A[i][k] * B[k][0];
     C[i][1] += A[i][k] * B[k][1];
     C[i][2] += A[i][k] * B[k][2];
     C[i][3] += A[i][k] * B[k][3];
  }

就是A[i][k]*B[k][0..3]

标量乘以向量。

.data
mat_A: .float 0.1, 0.2, 0.0, 0.1
       .float 0.2, 0.1, 0.3, 0.0
       .float 0.0, 0.3, 0.1, 0.5 
       .float 0.0, 0.6, 0.4, 0.1
mat_B: .float  4.92,  2.54, -0.63, -1.75
       .float  3.02, -1.51, -0.87,  1.35
       .float -4.29,  2.14,  0.71,  0.71
       .float -0.95,  0.48,  2.38, -0.95
mat_C: .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
 
format_result : .asciz "Matrix result is:\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n%5.3f %5.3f %5.3f %5.3f\n"
 
.text
naive_vec_matmul_4x4:
    /* r0 address of A
       r1 address of B
       r2 address of C
    */
push {r4, r5, r6, r7, r8, lr} /* 保存整型寄存器*/
    vpush {s16-s19}  //保存s16-s19
    vpush {s24-s27}

    mov r4, r2;//矩阵C地址给r4
    mov r5, #16;//矩阵C条目数量
    mov r6, #0;//矩阵C都赋值为0
    b .Lloop_init_test//初始化矩阵C
    .Lloop_init :
      str r6, [r4], +#4   /* *r4 ← r6 then r4 ← r4 + 4 */
    .Lloop_init_test:
      subs r5, r5, #1
      bge .Lloop_init
 
    /* 使用r4为i，r5为j，r6 为k*/
    mov r4, #0 /* r4 ← 0 ,k=0*/
    .Lloop_k:  /* 循环k,最外层 */
      cmp r4, #4  /* if r4 == 4 循环k结束 */
      beq .Lend_loop_k
      mov r5, #0  /* r5 ← 0,i从0开始 */
      .Lloop_i: /*循环j*/
       cmp r5, #4 /* if r5 == 4循环i结束 */
        beq .Lend_loop_i
        /* 计算C[i][0]地址并加载到s0 */
        /* C[i][0] 地址是 C + 4*(4 * i + 0) */
        add r7, r2, r5, LSL #4  
        vldmia r7, {s8-s15}   /* 矢量加载{s8,s9,s10,s11} ← {c[i][0], c[i][1], c[i][2], c[i][3]},c[i+1][0],...,c[i+1][3] */

        /* a[i][k]地址 */
        add r8, r4, r5,LSL #2
        add r8, r0, r8, LSL #2
        vldr s0, [r8]                  /* Load s0 ← a[i][k] */
        vldr s1, [r8,#16]                  /* Load s0 ← a[i+1][k] */

        /* b[k][0]地址 */
        add r8, r1,r4,LSL #4 /*r8 ← r1 + r4 << 4. r8 ← b + 4*(4*k)*/
        vldmia r8,{s16-s19}

        //vmul.f32 s24,s16,s0 /* {s24,s25,s26,s27} ← {s16,s17,s18,s19} * {s0}，就是a[i][k]*b[k][0…3]*/
        //vadd.f32 s8,s8,s24 /* {s8,s9,s10,s11} ← {s8,s9,s10,s11} + {s24,s25,s26,s27}


	vmla.f32 s8,s0,s16
	vmla.f32 s9,s0,s17
	vmla.f32 s10,s0,s18
	vmla.f32 s11,s0,s19

        vmla.f32 s12,s1,s16
	vmla.f32 s13,s1,s17
	vmla.f32 s14,s1,s18
	vmla.f32 s15,s1,s19
        vstmia r7,{s8-s15}      /* 保存s0 到 C[i][0]...C[i][3],c[i+1][0]...c[i+1][3] */

        add r5, r5, #2  /* r5 ← r5 + 1,j+2 */
        b .Lloop_i /*下个j迭代*/
       .Lend_loop_i: /* Here ends loop j */
       add r4, r4, #1 /* r4 ← r4 + 1，i++*/
       b .Lloop_k     /*下个i迭代 */
.Lend_loop_k: /* 循环k*/

    vpop {s24-s27}                
    vpop {s16-s19}
    pop {r4, r5, r6, r7, r8, lr}  /* Restore integer registers */
    bx lr /* Leave function */
 
.globl main
main:
    push {r4, r5, r6, lr}  /*保存整型寄存器*/
    vpush {d0-d1}          /* 保存浮点寄存器*/
 
    /* 准备调用函数naive_matmul_4x4 */
    ldr r0, addr_mat_A  /* r0 ← a */
    ldr r1, addr_mat_B  /* r1 ← b */
    ldr r2, addr_mat_C  /* r2 ← c */
    bl naive_vec_matmul_4x4
 
    /* 打印结果矩阵*/
    ldr r4, addr_mat_C  /* r4 ← c矩阵C地址 */
    vldr s0, [r4] /* s0 ← *r4. 第一个参数s0 ← c[0][0] */
    vcvt.f64.f32 d1, s0 /* 转成双精度d1 ← s0*/
    vmov r2, r3, d1      /* {r2,r3} ← d1 */
 
    mov r6, sp     /* 保存sp,r6 ← sp */
    mov r5, #1  /* 从1-15，共60个字节 */
    add r4, r4, #60 /* 矩阵c最后一个地址是c[3][3] */
    .Lloop:
        vldr s0, [r4] /* s0 ← *r4. 加载当前项*/
        vcvt.f64.f32 d1, s0 /*转成双精度d1 ← s0*/
        sub sp, sp, #8      /* 在堆栈预留一个双精度的空间*/
        vstr d1, [sp]       /* 将当前项保存到堆栈中 */
        sub r4, r4, #4      /* 矩阵地址减4个字节到前一项*/
        add r5, r5, #1      /* 待处理的矩阵项数目*/
        cmp r5, #16         /* 判断是会否结束*/
        bne .Lloop
 
    ldr r0, addr_format_result /* 格式化字符串*/
    bl printf /* 调用函数printf */
    mov sp, r6  /* 恢复sp寄存器*/
 
    mov r0, #0
    vpop {d0-d1}
    pop {r4, r5, r6, lr}
    bx lr
 
addr_mat_A : .word mat_A
addr_mat_B : .word mat_B
addr_mat_C : .word mat_C
addr_format_result : .word format_result

as -g -o matmul_best.o matmul_best.s -mfpu=neon

gcc -o matmul_best matmul_best.o

性能对比

测试工具：

.data
mat_A: .float 0.1, 0.2, 0.0, 0.1
       .float 0.2, 0.1, 0.3, 0.0
       .float 0.0, 0.3, 0.1, 0.5 
       .float 0.0, 0.6, 0.4, 0.1
mat_B: .float  4.92,  2.54, -0.63, -1.75
       .float  3.02, -1.51, -0.87,  1.35
       .float -4.29,  2.14,  0.71,  0.71
       .float -0.95,  0.48,  2.38, -0.95
mat_C: .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0


format_result : .asciz "Matrix result is:\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n"


.text
best_vectorial_matmul_4x4:
    /* r0 address of A
       r1 address of B
       r2 address of C
    */
push {r4, r5, r6, r7, r8, lr} /* 保存整型寄存器*/
    vpush {s16-s19}  //保存s16-s19
    vpush {s24-s27}

    mov r4, r2;//矩阵C地址给r4
    mov r5, #16;//矩阵C条目数量
    mov r6, #0;//矩阵C都赋值为0
    b .Lloop_init_test//初始化矩阵C
    .Lloop_init :
      str r6, [r4], +#4   /* *r4 ← r6 then r4 ← r4 + 4 */
    .Lloop_init_test:
      subs r5, r5, #1
      bge .Lloop_init
 
    /* 使用r4为i，r5为j，r6 为k*/
    mov r4, #0 /* r4 ← 0 ,k=0*/
    .Lloop_k:  /* 循环k,最外层 */
      cmp r4, #4  /* if r4 == 4 循环k结束 */
      beq .Lend_loop_k
      mov r5, #0  /* r5 ← 0,i从0开始 */
      .Lloop_i: /*循环j*/
       cmp r5, #4 /* if r5 == 4循环i结束 */
        beq .Lend_loop_i
        /* 计算C[i][0]地址并加载到s0 */
        /* C[i][0] 地址是 C + 4*(4 * i + 0) */
        add r7, r2, r5, LSL #4  
        vldmia r7, {s8-s15}   /* 矢量加载{s8,s9,s10,s11} ← {c[i][0], c[i][1], c[i][2], c[i][3]},c[i+1][0],...,c[i+1][3] */

        /* a[i][k]地址 */
        add r8, r4, r5,LSL #2
        add r8, r0, r8, LSL #2
        vldr s0, [r8]                  /* Load s0 ← a[i][k] */
        vldr s1, [r8,#16]                  /* Load s0 ← a[i+1][k] */

        /* b[k][0]地址 */
        add r8, r1,r4,LSL #4 /*r8 ← r1 + r4 << 4. r8 ← b + 4*(4*k)*/
        vldmia r8,{s16-s19}

        //vmul.f32 s24,s16,s0 /* {s24,s25,s26,s27} ← {s16,s17,s18,s19} * {s0}，就是a[i][k]*b[k][0…3]*/
        //vadd.f32 s8,s8,s24 /* {s8,s9,s10,s11} ← {s8,s9,s10,s11} + {s24,s25,s26,s27}


	vmla.f32 s8,s0,s16
	vmla.f32 s9,s0,s17
	vmla.f32 s10,s0,s18
	vmla.f32 s11,s0,s19

        vmla.f32 s12,s1,s16
	vmla.f32 s13,s1,s17
	vmla.f32 s14,s1,s18
	vmla.f32 s15,s1,s19
        vstmia r7,{s8-s15}      /* 保存s0 到 C[i][0]...C[i][3],c[i+1][0]...c[i+1][3] */

        add r5, r5, #2  /* r5 ← r5 + 1,j+2 */
        b .Lloop_i /*下个j迭代*/
       .Lend_loop_i: /* Here ends loop j */
       add r4, r4, #1 /* r4 ← r4 + 1，i++*/
       b .Lloop_k     /*下个i迭代 */
.Lend_loop_k: /* 循环k*/

    vpop {s24-s27}                
    vpop {s16-s19}
    pop {r4, r5, r6, r7, r8, lr}  /* Restore integer registers */
    bx lr /* Leave function */


.globl main
main:
    push {r4, r5, r6, lr}  /* Keep integer registers */

    ldr r0, addr_mat_A  /* r0 ← a */
    ldr r1, addr_mat_B  /* r1 ← b */
    ldr r2, addr_mat_C  /* r2 ← c */
    mov r4, #1
    mov r4, r4, LSL #21
    .Lmain_loop_test: 
      bl best_vectorial_matmul_4x4
      subs r4, r4, #1
      bne .Lmain_loop_test /* Should have been 'bge' */
    mov r0, #0
    pop {r4, r5, r6, lr}
    bx lr
addr_mat_A : .word mat_A
addr_mat_B : .word mat_B
addr_mat_C : .word mat_C
addr_format_result : .word format_result

as -g -o matmul_bench.o matmul_bench.s -mfpu=neon

gcc -o matmul_bench matmul_bench.o

通过#time ./matmul_bench得到时间

原始代码时间：

real 0m2.948s

向量优化后：

real 0m1.182s

寄存器并发后：

real 0m0.949s

访问优化后：

real 0m0.685s

相比原始的代码性能提升了4倍多。

64位矩阵乘

再来看下64位的机器上的表现，即ARMV8架构上的。

原始代码

.data
mat_A: .float 0.1, 0.2, 0.0, 0.1
       .float 0.2, 0.1, 0.3, 0.0
       .float 0.0, 0.3, 0.1, 0.5 
       .float 0.0, 0.6, 0.4, 0.1
mat_B: .float  4.92,  2.54, -0.63, -1.75
       .float  3.02, -1.51, -0.87,  1.35
       .float -4.29,  2.14,  0.71,  0.71
       .float -0.95,  0.48,  2.38, -0.95
mat_C: .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
 
format_result : .asciz "Matrix result is:\n%5.3f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.3f\n"
 
.text
naive_matmul_4x4:
    	/* x0 address of A
       	x1 address of B
       	x2 address of C */
	mov x4, x2;//矩阵C地址给r4
	mov x5, #16;//矩阵C条目数量
	mov x6, #0;//矩阵C都赋值为0
	b .Lloop_init_test//初始化矩阵C
	.Lloop_init :
	str x6, [x4], #+4   /* *x4 ← x6 then x4 ← x4 + 4 */
	.Lloop_init_test:
	subs x5, x5, #1
	bge .Lloop_init
 
	/* 使用x4为i，x5为j，x6 为k*/
	mov x4, #0 /* r4 ← 0 ,i=0*/
	.Lloop_i:  /* 循环i */
	cmp x4, #4  /* if x4 == 4 循环i结束 */
	beq .Lend_loop_i
	mov x5, #0  /* x5 ← 0,j从0开始 */
	.Lloop_j: /*循环j*/
	cmp x5, #4 /* if x5 == 4循环j结束 */
	beq .Lend_loop_j
	  /* 计算C[i][j]地址并加载到s0 */
	        /* C[i][j] 地址是 C + 4*(4 * i + j) */
          mov x7, x5               /* x7 ← x5,x7=j */
          adds x7, x7, x4, LSL #2  /* x7 ← x7 + (x4 << 2). 即x7 ← j + i * 4. We multiply i by the row size (4 elements) */
          adds x7, x2, x7, LSL #2  /* r7 ← r2 + (r7 << 2).即r7 ← C + 4*(j + i * 4)， 用原始大小乘以(j + i * 4) */
          ldr s0, [x7] /* s0 ← *r7 ，元素加载到s0,可设置断点*/
          mov x6, #0 /* r6 ← 0 ,k=0*/
          .Lloop_k :  /* 循环 k */
            cmp x6, #4 /* if x6 == 4循环 k 结束*/
            beq .Lend_loop_k
            /* 计算a[i][k]地址，并加载到s1*/
            /* a[i][k] 地址a + 4*(4 * i + k) */
            mov x8, x6  /* x8 ← x6. x8=k */
            adds x8, x8, x4, LSL #2  /* r8 ← r8 + (r4 << 2). 即r8 ← k + i * 4 */
            adds x8, x0, x8, LSL #2  /* r8 ← r0 + (r8 << 2). 即r8 ← a + 4*(k + i * 4) */
            ldr s1, [x8]    /* s1 ← *r8 ，即a[i][k]放到s1中，,可设置断点*/
            /* 计算b[k][j]地址，并加载到s2 */
            /* b[k][j] 地址b + 4*(4 * k + j) */
            mov x8, x5               /* r8 ← r5. 即r8 ← j */
            adds x8, x8, x6, LSL #2  /* r8 ← r8 + (r6 << 2). 即r8 ← j + k * 4 */
            adds x8, x1, x8, LSL #2  /* r8 ← r1 + (r8 << 2). 即r8 ← b + 4*(j + k * 4) */
            ldr s2, [x8]  /* s1 ← *r8 ，b[k][j]放到s2，,可设置断点*/
            fmul s3, s1, s2      /* s3 ← s1 * s2 */
            fadd s0, s0, s3      /* s0 ← s0 + s3 ，,可设置断点*/
 
            add x6, x6, #1           /* r6 ← r6 + 1,k++ */
           b .Lloop_k               /* 循环k下个迭代 */
          .Lend_loop_k: /* Here ends loop k */
        str s0, [x7]      /* 保存s0 回 C[i][j] ,可设置断点*/
        add x5, x5, #1  /* r5 ← r5 + 1,j++ */
        b .Lloop_j /*下个j迭代*/
       .Lend_loop_j: /* Here ends loop j */
       add x4, x4, #1 /* r4 ← r4 + 1，i++*/
       b .Lloop_i     /*下个i迭代 */
    .Lend_loop_i: /* Here ends loop i */
    ret
 
.globl _start
_start:
 
    /* 准备调用函数naive_matmul_4x4 */
    ldr x0, addr_mat_A  /* x0 ← a */
    ldr x1, addr_mat_B  /* x1 ← b */
    ldr x2, addr_mat_C  /* x2 ← c */
    bl naive_matmul_4x4
 
    /* 打印结果矩阵*/
    ldr x8, addr_mat_C  /* r4 ← c矩阵C地址 */
    ldr s0, [x8],#+4  /* s0 ← *r4. 第一个参数s0 ← c[0][0] */
    /* 转成双精度d1 ← s0*/
    fcvt d0,s0              /* 前八个参数从d0-d7 */
    ldr s1, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d1,s1
    ldr s2, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d2,s2
    ldr s3, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d3,s3
    ldr s4, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d4,s4
    ldr s5, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d5,s5
    ldr s6, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d6,s6
    ldr s7, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d7,s7
    mov x9, sp     /* 保存sp,r6 ← sp */
    ldr x11, addr_mat_C  /* x2 ← c */
    mov x10, #1  /* 从1-15，共60个字节 */
    add x11, x11, #32 /* 矩阵c的第九个参数地址4*8=32，给x11 */
    sub sp, sp, #80
    mov x12,sp //栈地址保存到x12中，这样复制过程sp可保持不变
    .Lloop:
        ldr  s8, [x11] /* 加载当前x11指向的矩阵条目*/
        fcvt d8, s8 /* 转成双精度*/
	str d8 ,[x12] /*加载到栈中*/
	add x12,x12,#8 /*栈往下移动，腾出下一个空位置*/
        add x11, x11, #4      /* 矩阵往下指一个*/
        add x10, x10, #1      /* 计数，共还需输出8个*/
        cmp x10, #9         /* 判断是会否结束，8个输出完毕*/
        bne .Lloop
 
    ldr x0, addr_format_result /* 格式化字符串*/
    bl printf /* 调用函数printf */
    mov sp, x9  /* 恢复sp寄存器*/

    mov x0, #0
    mov x8, 93
    svc 0

addr_mat_A : .dword mat_A
addr_mat_B : .dword mat_B
addr_mat_C : .dword mat_C
addr_format_result : .dword format_result

as -g -o matmul.o matmul.s

ld -o matmul matmul.o -lc -I /lib64/ld-linux-aarch64.so.1

这里的输出printf函数涉及到了大于8个参数的输出，printf第一个变量是通过x0传递，后面的可变参数前8个是通过d0-d7来传递，剩余通过栈来传递的。

64位向量优化

优化方法同32位架构中的代码，

由于A [i] [k]和A [i] [k + 1]是内存中的连续元素，因此访问A [i] [k]有多重加载。可以完全避免所有循环k并执行从A [i] [0]到A [i] [3]的4元素加载。

通过这种方法，可以一次删除4个操作，从而完全删除循环k。

此外，修改外层循环的次序，将K层移动到最外层，可以到的如下：

.data
mat_A: .float 0.1, 0.2, 0.0, 0.1
       .float 0.2, 0.1, 0.3, 0.0
       .float 0.0, 0.3, 0.1, 0.5 
       .float 0.0, 0.6, 0.4, 0.1
mat_B: .float  4.92,  2.54, -0.63, -1.75
       .float  3.02, -1.51, -0.87,  1.35
       .float -4.29,  2.14,  0.71,  0.71
       .float -0.95,  0.48,  2.38, -0.95
mat_C: .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
 
format_result : .asciz "Matrix result is:\n%5.3f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.3f\n"
 
.text
naive_matmul_4x4:
    	/* x0 address of A
       	x1 address of B
       	x2 address of C */
	mov x4, x2;//矩阵C地址给r4
	mov x5, #16;//矩阵C条目数量
	mov x6, #0;//矩阵C都赋值为0
	b .Lloop_init_test//初始化矩阵C
	.Lloop_init :
	str x6, [x4], #+4   /* *x4 ← x6 then x4 ← x4 + 4 */
	.Lloop_init_test:
	subs x5, x5, #1
	bge .Lloop_init
 
	/* 使用x4为i，x5为j，x6 为k*/
	mov x4, #0 /* r4 ← 0 ,i=0*/
	.Lloop_i:  /* 循环i */
	cmp x4, #4  /* if x4 == 4 循环i结束 */
	beq .Lend_loop_i
	mov x5, #0  /* x5 ← 0,j从0开始 */
	.Lloop_j: /*循环j*/
	cmp x5, #4 /* if x5 == 4循环j结束 */
	beq .Lend_loop_j
	  /* 计算C[i][j]地址并加载到s0 */
	        /* C[i][j] 地址是 C + 4*(4 * i + j) */
          mov x7, x5               /* x7 ← x5,x7=j */
          adds x7, x7, x4, LSL #2  /* x7 ← x7 + (x4 << 2). 即x7 ← j + i * 4. We multiply i by the row size (4 elements) */
          adds x7, x2, x7, LSL #2  /* r7 ← r2 + (r7 << 2).即r7 ← C + 4*(j + i * 4)， 用原始大小乘以(j + i * 4) */
          ldr q0,[x7];//矩阵c一行一次加载到q0
          mov x6, #0 /* r6 ← 0 ,k=0*/

            /* 计算a[i][0]地址，并加载到s1*/
            /* a[i][0] 地址a + 4*(4 * i + k) */
            mov x8, x6  /* x8 ← x6. x8=k */
            adds x8, x8, x4, LSL #2  /* r8 ← r8 + (r4 << 2). 即r8 ← k + i * 4 */
            adds x8, x0, x8, LSL #2  /* r8 ← r0 + (r8 << 2). 即r8 ← a + 4*(k + i * 4) */
            ldr q1, [x8]    
            /* 计算b[0][j]地址，并加载到s2 */
            /* b[0][j] 地址b + 4*(4 * k + j) */
            mov x8, x5               /* r8 ← r5. 即r8 ← j */
            adds x8, x8, x6, LSL #2  /* r8 ← r8 + (r6 << 2). 即r8 ← j + k * 4 */
            adds x8, x1, x8, LSL #2  /* r8 ← r1 + (r8 << 2). 即r8 ← b + 4*(j + k * 4) */
           ldr s16, [x8]       /* s16 ← *r8. 即s16 ← b[0][j] */
           ldr s17, [x8, #16] /* s17 ← *(r8 + 16). 即s17 ← b[1][j] */
           ldr s18, [x8, #32] /* s18 ← *(r8 + 32). 即s17 ← b[2][j] */
           ldr s19, [x8, #48] /* s19 ← *(r8 + 48). 即s17 ← b[3][j] */
	   sub sp, sp,#64 //准备栈，先保存到栈中，然后再倒入到q2中。
	   mov x9,sp
	    st4 {v16.s,v17.s,v18.s,v19.s}[0],[x9]//导入到栈中
            ldr q2,[x9]// 导入到q2中。
            fmul v3.4s, v1.4s, v2.4s     /* 向量矩阵乘，行乘以列*/
	    str q3,[x9] //结果倒回到栈中，因为还要将q3中的3个单精度数进行相加。
            ldr s4,[x9,#+4]//第二个单精度数，第一个就是s3
	    fadd s0,s3,s4
            ldr s4,[x9,#+8] //第三个单精度数
	    fadd s0, s0,s4
            ldr s4,[x9,#+12] //第四个单精度数
	    fadd s0, s0,s4
	     
	    add sp,sp,#64//恢复栈

        str s0, [x7]      /* 保存s0 回 C[i][j] ,可设置断点*/
        add x5, x5, #1  /* r5 ← r5 + 1,j++ */
        b .Lloop_j /*下个j迭代*/
       .Lend_loop_j: /* Here ends loop j */
       add x4, x4, #1 /* r4 ← r4 + 1，i++*/
       b .Lloop_i     /*下个i迭代 */
    .Lend_loop_i: /* Here ends loop i */
    ret
 
.globl _start
_start:
 
    /* 准备调用函数naive_matmul_4x4 */
    ldr x0, addr_mat_A  /* x0 ← a */
    ldr x1, addr_mat_B  /* x1 ← b */
    ldr x2, addr_mat_C  /* x2 ← c */
    bl naive_matmul_4x4
 
    /* 打印结果矩阵*/
    ldr x8, addr_mat_C  /* r4 ← c矩阵C地址 */
    ldr s0, [x8],#+4  /* s0 ← *r4. 第一个参数s0 ← c[0][0] */
    /* 转成双精度d1 ← s0*/
    fcvt d0,s0              /* 前八个参数从d0-d7 */
    ldr s1, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d1,s1
    ldr s2, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d2,s2
    ldr s3, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d3,s3
    ldr s4, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d4,s4
    ldr s5, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d5,s5
    ldr s6, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d6,s6
    ldr s7, [x8],#+4       /* 将当前项保存到堆栈中 */
    fcvt d7,s7
    mov x9, sp     /* 保存sp,r6 ← sp */
    ldr x11, addr_mat_C  /* x2 ← c */
    mov x10, #1  /* 从1-15，共60个字节 */
    add x11, x11, #32 /* 矩阵c的第九个参数地址4*8=32，给x11 */
    sub sp, sp, #80
    mov x12,sp //栈地址保存到x12中，这样复制过程sp可保持不变
    .Lloop:
        ldr  s8, [x11] /* 加载当前x11指向的矩阵条目*/
        fcvt d8, s8 /* 转成双精度*/
	str d8 ,[x12] /*加载到栈中*/
	add x12,x12,#8 /*栈往下移动，腾出下一个空位置*/
        add x11, x11, #4      /* 矩阵往下指一个*/
        add x10, x10, #1      /* 计数，共还需输出8个*/
        cmp x10, #9         /* 判断是会否结束，8个输出完毕*/
        bne .Lloop
 
    ldr x0, addr_format_result /* 格式化字符串*/
    bl printf /* 调用函数printf */
    mov sp, x9  /* 恢复sp寄存器*/

    mov x0, #0
    mov x8, 93
    svc 0

addr_mat_A : .dword mat_A
addr_mat_B : .dword mat_B
addr_mat_C : .dword mat_C
addr_format_result : .dword format_result

as -g -o matmulo1.o matmulo1.s

ld -o matmulo1 matmulo1.o -lc -I /lib64/ld-linux-aarch64.so.1

tips：其中fmla是a=a+b×c

64位性能对比

.data
mat_A: .float 0.1, 0.2, 0.0, 0.1
       .float 0.2, 0.1, 0.3, 0.0
       .float 0.0, 0.3, 0.1, 0.5 
       .float 0.0, 0.6, 0.4, 0.1
mat_B: .float  4.92,  2.54, -0.63, -1.75
       .float  3.02, -1.51, -0.87,  1.35
       .float -4.29,  2.14,  0.71,  0.71
       .float -0.95,  0.48,  2.38, -0.95
mat_C: .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0
       .float 0.0, 0.0, 0.0, 0.0

format_result : .asciz "Matrix result is:\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n%5.2f %5.2f %5.2f %5.2f\n"

.text
naive_matmul_4x4:
    	/* x0 address of A
       	x1 address of B
       	x2 address of C */
	mov x4, x2;//矩阵C地址给r4
	mov x5, #16;//矩阵C条目数量
	mov x6, #0;//矩阵C都赋值为0
	b .Lloop_init_test//初始化矩阵C
	.Lloop_init :
	str x6, [x4], #+4   /* *x4 ← x6 then x4 ← x4 + 4 */
	.Lloop_init_test:
	subs x5, x5, #1
	bge .Lloop_init
 
	/* 使用x4为i，x5为j，x6 为k*/
	mov x4, #0 /* r4 ← 0 ,i=0*/
	.Lloop_i:  /* 循环i */
	cmp x4, #4  /* if x4 == 4 循环i结束 */
	beq .Lend_loop_i
	mov x5, #0  /* x5 ← 0,j从0开始 */
	.Lloop_j: /*循环j*/
	cmp x5, #4 /* if x5 == 4循环j结束 */
	beq .Lend_loop_j
	  /* 计算C[i][j]地址并加载到s0 */
	        /* C[i][j] 地址是 C + 4*(4 * i + j) */
          mov x7, x5               /* x7 ← x5,x7=j */
          adds x7, x7, x4, LSL #2  /* x7 ← x7 + (x4 << 2). 即x7 ← j + i * 4. We multiply i by the row size (4 elements) */
          adds x7, x2, x7, LSL #2  /* r7 ← r2 + (r7 << 2).即r7 ← C + 4*(j + i * 4)， 用原始大小乘以(j + i * 4) */
          ldr s0, [x7] /* s0 ← *r7 ，元素加载到s0,可设置断点*/
          mov x6, #0 /* r6 ← 0 ,k=0*/
          .Lloop_k :  /* 循环 k */
            cmp x6, #4 /* if x6 == 4循环 k 结束*/
            beq .Lend_loop_k
            /* 计算a[i][k]地址，并加载到s1*/
            /* a[i][k] 地址a + 4*(4 * i + k) */
            mov x8, x6  /* x8 ← x6. x8=k */
            adds x8, x8, x4, LSL #2  /* r8 ← r8 + (r4 << 2). 即r8 ← k + i * 4 */
            adds x8, x0, x8, LSL #2  /* r8 ← r0 + (r8 << 2). 即r8 ← a + 4*(k + i * 4) */
            ldr s1, [x8]    /* s1 ← *r8 ，即a[i][k]放到s1中，,可设置断点*/
            /* 计算b[k][j]地址，并加载到s2 */
            /* b[k][j] 地址b + 4*(4 * k + j) */
            mov x8, x5               /* r8 ← r5. 即r8 ← j */
            adds x8, x8, x6, LSL #2  /* r8 ← r8 + (r6 << 2). 即r8 ← j + k * 4 */
            adds x8, x1, x8, LSL #2  /* r8 ← r1 + (r8 << 2). 即r8 ← b + 4*(j + k * 4) */
            ldr s2, [x8]  /* s1 ← *r8 ，b[k][j]放到s2，,可设置断点*/
            fmul s3, s1, s2      /* s3 ← s1 * s2 */
            fadd s0, s0, s3      /* s0 ← s0 + s3 ，,可设置断点*/
 
            add x6, x6, #1           /* r6 ← r6 + 1,k++ */
           b .Lloop_k               /* 循环k下个迭代 */
          .Lend_loop_k: /* Here ends loop k */
        str s0, [x7]      /* 保存s0 回 C[i][j] ,可设置断点*/
        add x5, x5, #1  /* r5 ← r5 + 1,j++ */
        b .Lloop_j /*下个j迭代*/
       .Lend_loop_j: /* Here ends loop j */
       add x4, x4, #1 /* r4 ← r4 + 1，i++*/
       b .Lloop_i     /*下个i迭代 */
    .Lend_loop_i: /* Here ends loop i */
    ret
.globl main
main:

    /* 准备调用函数naive_matmul_4x4 */
    ldr x0, addr_mat_A  /* x0 ← a */
    ldr x1, addr_mat_B  /* x1 ← b */
    ldr x2, addr_mat_C  /* x2 ← c */

    mov x4, #1
    mov x4, x4, LSL #25
    sub sp,sp,#16
    .Lmain_loop_test: 
      str x4,[sp]
      bl naive_matmul_4x4
      ldr x4,[sp]
      sub x4, x4, #1
     cmp x4,#1
     bge .Lmain_loop_test /* Should have been 'bge' */
    add sp,sp,#16
    mov x0, #0
    mov x8, 93
    svc 0

addr_mat_A : .dword mat_A
addr_mat_B : .dword mat_B
addr_mat_C : .dword mat_C
addr_format_result : .dword format_result

as -g -o matmul_bench.o matmul_bench.s

gcc -o matmul_bench matmul_bench.o

通过#time ./matmul_bench得到时间

real 0m5.913s

使用优化过的需要时间：

real 0m4.088s

性能提升30%左右。

badman250

关注

0
点赞
踩
4

收藏

觉得还不错? 一键收藏
0
评论
第26部分- Linux ARM汇编浮点和向量

第26部分- Linux ARM汇编浮点和向量二进制浮点数是一个实数的近似表示，由三个部分组成：符号，尾数和指数。1.01110是20+ 2-2+ 2-3+ 2-4= 1.43750(10)VFPv2支持两个IEEE 754数字：Binary32和Binary64，通常以其C类型分别为float和double或单精度和双精度。在单精度浮点中，尾数为23位（归一化数字为整数1的+1），指数为8位（因此指数范围为-126到127）。在双精度浮点中，尾数为52位（+1），指数...
复制链接

扫一扫