矩阵分块乘法

矩阵乘法

普通矩阵乘法（ijk）

void matrix_multply(float *a, float *b, float *c,  int n)
{
    for(int i = 0; i < n; i++)
    {
        for(int j = 0; j < n; j++)
        {
            float sum = c[i*N + j] ;
            for(int k = 0; k < n; k++)
            {
                sum += a[i*N + k] * b[k*N + j];
            }
            c[i*N + j] = sum; 
        }
    }        
}

快速矩阵乘法（ikj）
跳转更少，cache命中率更高

void matrix_multply_quick(float *a, float *b, float *c,  int n)
{
    float s;
    for(int i = 0; i < n; i++){
        for(int k = 0; k < n; k++){
            s = a[i*N + k];
            for(int j = 0; j < n; j++){
                c[i*N + j] += s*b[k*N + j];
            }
        }
    }
}

https://zhuanlan.zhihu.com/p/146250334

矩阵分块乘法

void dgemm_block (float* a, float* b, float* c, int n)
{
    for ( int sj = 0; sj < n; sj += M )

        for ( int si = 0; si < n; si += M )

            for ( int sk = 0; sk < n; sk += M )

                matrix_multply(a+si*n+sk, b+sk*n+sj,  c+si*n+sj, M);

}

https://zhuanlan.zhihu.com/p/129392257?utm_source=wechat_session&utm_medium=social&utm_oi=807246749705080832&utm_campaign=shareopn

• 对A的引用具有良好的空间局部性，因为每个行片段的迭代的步长为1；与此同时，还有很好的时间局部性，因为整个行片段的棉条被连续引用了bsize次。

• 对B的引用具有良好的时间局部性，因为整个bsize×bsize块被连续访问n次。

• 对C的引用具有良好的空间局部性，因为行片段的每个元素都是被连续写入的。
  请注意，对C的引用没有良好的时间局部性，因为每个行片段只能访问一次。

cache

http://cenalulu.github.io/linux/all-about-cpu-cache/#:~:text=%E4%BB%80%E4%B9%88%E6%98%AFCache%20Line%20Cache%20Line%E5%8F%AF%E4%BB%A5%E7%AE%80%E5%8D%95%E7%9A%84%E7%90%86%E8%A7%A3%E4%B8%BACPU,Cache%E4%B8%AD%E7%9A%84%E6%9C%80%E5%B0%8F%E7%BC%93%E5%AD%98%E5%8D%95%E4%BD%8D%E3%80%82%E7%9B%AE%E5%89%8D%E4%B8%BB%E6%B5%81%E7%9A%84CPU%20Cache%E7%9A%84Cache%20Line%E5%A4%A7%E5%B0%8F%E9%83%BD%E6%98%AF64Bytes%E3%80%82%E5%81%87%E8%AE%BE%E6%88%91%E4%BB%AC%E6%9C%89%E4%B8%80%E4%B8%AA512%E5%AD%97%E8%8A%82%E7%9A%84%E4%B8%80%E7%BA%A7%E7%BC%93%E5%AD%98%EF%BC%8C%E9%82%A3%E4%B9%88%E6%8C%89%E7%85%A764B%E7%9A%84%E7%BC%93%E5%AD%98%E5%8D%95%E4%BD%8D%E5%A4%A7%E5%B0%8F%E6%9D%A5%E7%AE%97%EF%BC%8C%E8%BF%99%E4%B8%AA%E4%B8%80%E7%BA%A7%E7%BC%93%E5%AD%98%E6%89%80%E8%83%BD%E5%AD%98%E6%94%BE%E7%9A%84%E7%BC%93%E5%AD%98%E4%B8%AA%E6%95%B0%E5%B0%B1%E6%98%AF512%2F64%20%3D%208