为什么转置一个512x512的矩阵，会比513x513的矩阵慢很多？

谨以此文，纪念刚退休的Professor Sibert以及Professor Goel。你们尽管已年过70，却还仍然坚持在教导学生，实在令人钦佩。我今天所拥有的编程知识，经验，技巧，很大一部分是从你们那儿学来的。谢谢你们。

问题出处：http://stackoverflow.com/questions/11413855/why-is-transposing-a-matrix-of-512x512-much-slower-than-transposing-a-matrix-of

事情的起因是这样的，先看下面这段代码：

#define SAMPLES 1000 #define MATSIZE 512

#include<time.h> #include<iostream> int mat[MATSIZE][MATSIZE];

void transpose() { for(int i =0; i < MATSIZE ; i++) for(int j =0; j < MATSIZE ; j++) { int aux = mat[i][j]; mat[i][j]= mat[j][i]; mat[j][i]= aux; } }

int main() { //initialize matrix for(int i =0; i < MATSIZE ; i++) for(int j =0; j < MATSIZE ; j++) mat[i][j]= i+j;

int t = clock(); for(int i =0; i < SAMPLES ; i++) transpose(); int elapsed = clock()- t;

std::cout <<"Average for a matrix of "<< MATSIZE <<": "<< elapsed / SAMPLES;

很普通的一个求矩阵转置的程序。但是，当MATSIZE取512和513的时候，出现了非常有意思的结果：

512 平均 2.19ms

513 平均 0.57ms

很让人惊讶吧，513竟然比512快。更进一步的研究发现，size=512的时候，运算速度会比同数量级的其它数字慢很多很多。这是怎么一回事呢？

stackoverflow上大牛给的解答非常正确，但是这次，我不想做翻译了。我从Professor Sibert那里，从Professor Goel那里，学到的知识，足够帮我解决这个问题了，我不是一个人。下面是我的解答：

很容易就联想到，造成这个问题的原因是CPU cache，我们有很多种方式来存储cache，具体可以参考这里。

原作者没有给出他的CPU型号，但是如今的pc几乎都是采用的set associative的cache结构，下面我用2-way set associate来做例子，讲解一下cache的工作原理。

（图片取自Professor Sibert的讲义，这可是纯ascii画的哦= =）

一个内存地址，可以划分为block，tag，word，byte 4个部分。10bits的block，对应了1024个cache set，内存地址的block固定了，就必须存储在相应的set里面，这样可以把查询cache的事件从O(n)缩短为O(1)。

举个例子，block是1023（1111111111），你的数据就放在第1023个set里面。可能有人会觉得奇怪，为什么block不是取的最前面的10bits，这当然是有道理的，通常在内存里数据都是连续存放的，就是说，同一段程序用的数据，他们前10位几乎都是一样的，如果用前10位来定位block，那么collision的发生率非常高，cache效率非常低下，所以才选了后面的10位来定位block。

当然，每个set里面有多条记录，2-way是2条，你得遍历这两条记录，比较前面50位的tag，如果tag一样，并且Valid bit(V)=1，那么恭喜你，你的数据在cache里面，接着就可以通过word和byte来取数据了。

如果遍历完这两条记录，还是没有找到tag的话，那么很遗憾，你的数据不在cache里，得从内存里读。从内存里获取相应的数据，然后把它存到对应的cache set里，如果set里有空位的话最好，如果没有的话，用LRU来替换。因为一个set里只有2条数据，所以实现LRU仅仅需要一个额外bit就可以了，非常高效。

好了，背景知识介绍的差不多了，让我们回到这个问题上来。为什么512大小的矩阵，会比其它数字慢那么多？

让我们来计算一下，512x512的int矩阵，在内存里是连续存放的。每个cache line是16bytes，对应4个int，所以一个n阶矩阵的row可以填充n/4个cache set。假设第一个数据a[0][0]正好对应cache set 0，那么其中每一个数据a[i][j]对应的cache set是(512i+j)/4%1024=(128i+j/4)%1024。可以看到，前面的系数正好可以整除。很不巧的是，在进行矩阵转置的运算时，在第2个for循环中，我们需要依次访问每一个row中对应i的值。这样会造成下面的结果：假设i=0，set(a[0][0])=0, set(a[1][0])=128, set(a[2][0])=256...set(a[7][0])=896,set(a[0][0])=0,后面开始重复了，到a[15][0]的时候刚好填完整个cache的所有128整数倍的set，当读取a[16][0]的时候，将会发生replace，把a[0][0]从cache里移除。这样，当源程序的i=1时，将完全重复i=0的计算过程，每次取数据都需要先从memory读到cache中来，cache的作用完全没有体现。

而当size=513的时候，事情就不一样了，mat[i][j]对应的cache set是(513i+j)/4%1024，前面的系数除不尽了，每递增4次结果会比size=512时偏差1。例如：set(a[0][0])=0, set(a[1][0])=128, set(a[2][0])=256，set(a[3][0])=384, set(a[4][0])=513...这样就很微妙的把cache set给错开了。a[16][0]不在第0行而是第4行，不会覆盖之前的数据。即使将全部的a[0-15][i]都读入cache，也不会发生碰撞。之后，由于一个cache有4个word，a[0-15][i+1]，a[0-15][i+2]，a[0-15][i+3]也同时被读进cache里了，所以计算i+1,i+2,i+3时，仅仅需要读对应行的数据就可以了，同一行的数据都是连续的，所以碰撞率很低。这个计算过程很好的利用了cache，如果不考虑其他因素（实际上，这个已经是影响运行时间的最大因素了），理论上我们可以节省75%的运行时间，可以看到，这个理论预测是和提问者给的数据相符合的。

总之，当你的data size是128的整数倍的时候，得特别小心，搞不好cache collision就把你的程序给拖慢了呢

Update 1: 原代码有逻辑错误，这点大家都不要吐槽了，代码不是我写的= =

Update 2：帅哥问我，为什么可以加速这么多。这个循环包括4次读cache的操作，2次写cache的操作，以及0-2次replace操作。每次replace操作会有一次memory read，有可能会有memory write（假设它是write back）。前面的读写cache时间和读写内存相比，几乎可以忽略，对效率产生显著影响的是后面的内存读写。如果cache的hit率高了，那么内存读写的次数就少了，程序运行时间是会产生很大影响的

Update 3：当然，具体效果还视乎CPU架构而定，我自己试验的只有节省25%左右时间

Update 4: 有人提出了用划分矩阵（把大矩阵分成若干个小矩阵分别计算）的方法来求转置。划分矩阵可以解决类似的问题（譬如说求两个矩阵乘积），但是对解决这个问题没有任何帮助。因为求转置的时候，每个数据只用到了一次，没有重复访问；即便划分成更小的矩阵，在cache里面的位置也没有发生改变。

Update 5: 据说，Professor Goel只是因病休息几个学期，没有退休。。。（原来你还要回来教课！！！）

为什么转置一个512x512的矩阵，会比513x513的矩阵慢很多？相关推荐

矩阵的基础知识回顾:矩阵乘法，矩阵的逆，伴随矩阵，矩阵的转置，行列式，相似矩阵，实对称矩阵
Agenda 1. 矩阵matrix 1.1 矩阵运算matrix operations 1.1.1 矩阵乘法matrix multiplication 1.1.1.1 简化矩阵乘法(facilita ...
【MATLAB实验】MATLAB矩阵与数组及改变矩阵形状（rot90函数逆时针旋转、矩阵转置）
目录 Matlab矩阵矩阵除法矩阵乘方数组的乘和除数组的乘方数据的输出格式常用函数的应用矩阵的建立冒号表达式结构矩阵和单元矩阵结构矩阵: 单元矩阵: 矩阵元素的引用方式利用冒号表 ...
数组：矩阵快速转置矩阵相加三元组顺序表/三元矩阵随机生成稀疏矩阵压缩矩阵【C语言，数据结构】（内含源代码）
目录题目: 题目分析: 概要设计: 二维矩阵数据结构: 三元数组\三元顺序表顺序表结构: 详细设计: 三元矩阵相加: 三元矩阵快速转置: 调试分析: 用户手册: 测试结果: 源代码: 主程序: 头文 ...
Java黑皮书课后题第8章：**8.14（探讨矩阵）编写程序，提示用户输入一个方阵的长度，随机地在矩阵中填入0和1，打印这个矩阵，然后找出整行、整列或者对角线都是1或0的行、列和对角线
**8.14(探讨矩阵)编写程序,提示用户输入一个方阵的长度,随机地在矩阵中填入0和1,打印这个矩阵,然后找出整行.整列或者对角线都是1或0的行.列和对角线题目题目描述与运行示例破题代码题目 ...
一个N*M的矩阵，找出这个矩阵中所有元素的和不小于K的面积最小的子矩阵
题目描述: 一个N*M的矩阵,找出这个矩阵中所有元素的和不小于K的面积最小的子矩阵(矩阵中元素个数为矩阵面积) 输入: 每个案例第一行三个正整数N,M<=100,表示矩阵大小,和一个整数K 接下 ...
螺旋矩阵一个n行n列的螺旅矩阵可由如下方法生成：从矩阵的左上角（第1行第 1列）出发，初始时向右移动；如果前方是未曾经过的格子，则继续前进，否则右转；重复上述操作直至经过矩阵中所有格子。
螺旋矩阵一个n行n列的螺旅矩阵可由如下方法生成: 从矩阵的左上角(第1行第 1列)出发,初始时向右移动:如果前方是未曾经过的格子,则继续前进,否则右转:重复上述操作直至经过矩阵中所有格子.根据经过顺 ...
python如何求矩阵逆运算_python简单实现矩阵的乘，加，转置和逆运算示例
本文实例讲述了python简单实现矩阵的乘,加,转置和逆运算.分享给大家供大家参考,具体如下: 使用python完成矩阵的乘,加,转置和逆: # -*- coding:utf-8 -*- #矩阵的乘, ...
MATLAB笔记5：矩阵的转置、求逆、旋转、翻转；矩阵的行列式、秩、迹；矩阵的特征值、特征向量
矩阵的转置.求逆.旋转.翻转 inv(A):求矩阵A的逆矩阵: 转置:A.'为矩阵A的转置,A'为矩阵A的共轭转置: rot90(A,k):将矩阵A逆时针方向旋转90°的k倍,k为1时可省略: fli ...
用Python编程创建一个10*10的ndarray对象，且矩阵边界全部为1，里面全部为0
一.问题描述用Python编程创建一个10*10的ndarray对象,且矩阵边界全部为1,里面全部为0 二.问题解答代码: import numpy as np a= np.zeros((10,1 ...

为什么转置一个512x512的矩阵，会比513x513的矩阵慢很多？

为什么转置一个512x512的矩阵，会比513x513的矩阵慢很多？相关推荐

最新文章

热门文章