HNU-2024计算机系统-小班讨论5

写在前面:
本次小班的选题依旧开放性很强，需要通过阅读老师发的pdf文档来完成相应的程序设计，从而实现矩阵乘法的高性能计算，通过本次讨论课的学习，同学们可以对程序性能优化有一个更加清楚的认识和理解。

本文所有题目由全组同学共同完成，感谢每一位同学的支持与付出。

一、选题

选题一
请阅读所附《高性能计算:矩阵乘法》pdf文档，设计一段不小于 30*30矩阵的乘法运算 C代码，请尝试使用你们的代码进行该文档中提及的优化方式，根据学习和实践所得，讨论为什么会确实优化了性能(并尝试分析各种优化后的 Cache 的命中率/不命中率)并将你们的讨论成果与全班分享。
Tips：这是一道具有综合应用属性的讨论选题，请一定要展示你们的优化过程和结果，并对所得结果进行尽可能的分析“为什么就加快了？”

二、分析与解答

1.优化方法一：改变计算顺序

2.优化方法二：采用分块矩阵

3.优化方法三：使用SIMD指令集进行计算优化

4.SIMD指令集的代码：

#include<immintrin.h>
#include <emmintrin.h> 
#include<stdio.h> 
#include <time.h>
#define N 32
double A[N][N];
double B[N][N];
double C[N][N];
int judge=1;
void Gemm_4x4Kernel(int k, double *A, int ldA, double *B, int ldB, double *C, int ldc) {
    // 创建一个向量寄存器，用于储存C的一行
    __m256d gamma_0123_0 = _mm256_loadu_pd(&C[0 + 0*ldc]);
    __m256d gamma_0123_1 = _mm256_loadu_pd(&C[0 + 1*ldc]);
    __m256d gamma_0123_2 = _mm256_loadu_pd(&C[0 + 2*ldc]);
    __m256d gamma_0123_3 = _mm256_loadu_pd(&C[0 + 3*ldc]);
    for (int p = 0; p < k; p++) {
        // 创建一个向量寄存器，用于储存C的一行
        __m256d alpha_0123_p = _mm256_broadcast_sd(&A[p + 0*ldA]);
        __m256d beta_p_j = _mm256_loadu_pd(&B[0 + p*ldB]);
        // 更新C中的值
        gamma_0123_0 = _mm256_fmadd_pd(alpha_0123_p, beta_p_j, gamma_0123_0);
        alpha_0123_p = _mm256_broadcast_sd(&A[p + 1*ldA]);
        gamma_0123_1 = _mm256_fmadd_pd(alpha_0123_p, beta_p_j, gamma_0123_1);
        alpha_0123_p = _mm256_broadcast_sd(&A[p + 2*ldA]);
        gamma_0123_2 = _mm256_fmadd_pd(alpha_0123_p, beta_p_j, gamma_0123_2);
        alpha_0123_p = _mm256_broadcast_sd(&A[p + 3*ldA]);
        gamma_0123_3 = _mm256_fmadd_pd(alpha_0123_p, beta_p_j, gamma_0123_3);
    }

    // 将寄存器的值写回数组
    _mm256_storeu_pd(&C[0 + 0*ldc], gamma_0123_0);
    _mm256_storeu_pd(&C[0 + 1*ldc], gamma_0123_1);
    _mm256_storeu_pd(&C[0 + 2*ldc], gamma_0123_2);
    _mm256_storeu_pd(&C[0 + 3*ldc], gamma_0123_3);
}
void matmul(double A[N][N], double B[N][N], double C[N][N]) {
    for (int i = 0; i < N; i += 4) {
        for (int j = 0; j < N; j += 4) {
            for (int k = 0; k < N; k += 4) {
                Gemm_4x4Kernel(4, &A[i][k], N, &B[k][j], N, &C[i][j], N);
            }
        }
    }
}
int main(){
    clock_t start, end;
	for(int i=0;i<32;i++){
		for(int j=0;j<32;j++){
			scanf("%lf",&A[i][j]);
		}
	}
	for(int i=0;i<32;i++){
		for(int j=0;j<32;j++){
			scanf("%lf",&B[i][j]);
		}
	}
    start = clock();
	matmul(A,B,C);
    end = clock();
	for(int i=0;i<32;i++){
		for(int j=0;j<32;j++){
			printf("%lf ",C[i][j]);
		}
        printf("\n");
	}
    printf("\nTime（ms） -> %f\n\n", (float)(end - start) * 1000 / CLOCKS_PER_SEC);
	return 0;
}