关于编译器优化方面的知识

mathe

对于循环语句，数据依赖关系可以根据它同循环体的关系，又分成循环内和跨循环体的依赖关系，比如

1. for(i=0;i<n;i++){
   
2.     X = A;                  (1)
   
3.     Y = X+2;                (2)
   
4.     A= 3*Y-1;               (3)
   
5. }

复制代码

对于上面的代码，语句(2)真依赖于语句(1),而这种依赖是在一次循环迭代内部的，这种在一次循环迭代内部的依赖我们称为循环内的依赖关系(intra-loop dependence).同样语句(1)真依赖于语句(3),但是这种依赖关系是跨越迭代的边界的，比如这里每次迭代的语句(1)用到了上次迭代时语句(3)产生的结果，对于这种跨越迭代边界的数据依赖关系，我们称为跨循环体的迭代关系(loop-carried dependence)

对于跨循环体的依赖关系，通常如果所有的读写都发生在某个数组上时，我们可以准确计算出这种数据依赖关系到底跨越了多少次迭代。比如

1. for(i=0;i<n;i++){
   
2.      for(j=0;j<n;j++){
   
3.            …=a[i][j]+…;     (1)
   
4.            a[i+1][j+2]=…..  (2)
   
5.      }
   
6. }

复制代码

我们可以知道i=0,j=0时语句(2)定义的内容将被i=1,j=2时的语句(1)所使用
更一般的，i=X,j=Y时语句(2)定义的内容将被i=X+1,j=Y+2时的语句(1)所使用。
我们计算(X+1,Y+2)-(X,Y)=(1,2).我们发现这个差值是一个常向量，同具体定义和使用这个内存空间的语句没有关系。这个常向量我们称之为这个数据依赖关系的依赖向量(Dependence Vector).
计算依赖向量时需要注意的是它们并不是对于的数组下标的差值。

而么模变换就是以依赖向量作为分析的基础。

mathe

现在开始考虑前面的问题，对于一个给定的循环，采用下列的循环优化如：
循环交换(Loop Interchange),循环反向(Loop Reversal),循环扭曲(Loop Skewing)等变化以及它们的组合
是否合法。
这个就是么模变换要讨论的问题。
么模变换首先假设循环体中的所有数据依赖关系都可以表示成依赖向量。这是一个很强的前提条件，通常，这个
要求循环中牵涉到的数据都是数组，不然，数据依赖关系就很难表达成依赖向量了。所以通常适合么模变换的代码
都是科学计算方面的代码。此外，编译器还要事先对代码做一些简单的变换，比如下面的代码：

1. for(i=0;i<n;i++){
   
2.    for(j=0;j<n;j++){
   
3.           x= A[i][j];
   
4.           B[i][j]+=x;
   
5.    }
   
6. }

复制代码

显然可以看出，对于变量x,所有的循环迭代之间都有关于x的输出依赖关系，所以这种依赖关系就无法表达成依赖向量，
但是我们可以通过消去标量x，将代码变化成

1. for(i=0;i<n;i++){
   
2.    for(j=0;j<n;j++){
   
3.           B[i][j]+= A[i][j];
   
4.    }
   
5. }

复制代码

那么就没有这个问题了，现在唯一的数据依赖就是依赖向量(0,0).
同样对于一些无法消除的关于标量的数据依赖关系，也可以通过将标量转化成数组来达到我们的目的，但是这种转化通常情况代价
比较大，很少会真正采用，比如上面的代码也可以转化为一下等价的代码：

1. for(i=0;i<n;i++){
   
2.    for(j=0;j<n;j++){
   
3.           X[i][j]= A[i][j];
   
4.           B[i][j]+=X[i][j];
   
5.    }
   
6. }

复制代码

但是显然这种转化不是很好。

此外还有一种特殊情况我们可以利用，这中情况在实际应用中会经常遇上，比如代码

1. for(i=0;i<n;i++){
   
2.    ...
   
3.    err+=(A[i]-B[i])*(A[i]-B[i]);
   
4.    if(...){
   
5.           err -= A[i]*B[i];
   
6.    }
   
7. }

复制代码

对于上面这个代码,显然关于标量err,同样所有的循环迭代之间都有输出依赖关系，这种依赖关系无法表达成依赖向量。
但是我们知道，如果上面的代码中，所有对于err的操作都是+=和-=某个同err无关的数据，由于我们知道+和-法操作满足结合率，
也就是我们可以任意改变上面所有累加累减的顺序而不会改变代码的正确性，所以实际上，这种依赖关系在我们的循环变换中
也可以忽略掉。所以对于这样的代码，我们依旧可以采用么模变换的判断方法。
想上面这种由于满足结合率导致我们可以忽略其数据依赖关系的操作符，我们称为Reduce操作符。典型的Reduce操作符有+,-,*,/(浮点除法，
整数除法由于截断，不符合要求),移位，逻辑运算，位运算等。
而最常见的就是+,-运算。通常的编译器都应该能够支持将整数的+,-运算作为Reduce操作符。而关于普通整数乘法由于实际上很难遇上所以一般
不会采用。而对于浮点数的乘除运算，通常编译器由于考虑到会对计算精度产生影响而不采用（对于部分编译器打开特殊优化选项则会采用）。

有了以上分析，我们可以大概知道对于那些循环可以采用么模变换的方法来判断一个循环变换是否合法了。
首先，这个循环中所有数据依赖关系都可以写成一批依赖向量，
我们可以注意到，所有这些依赖向量都必然是非负的(0或正的),这里正负的定义采用我们前面定义的方法（也就是第一项非零数据为正的称为正向量）。
我们继续查看前面用到的一个关于循环变换的图片

可以看出，每个变换都被一个矩阵所标识，实际上，由于循环变换而引起循环迭代空间上发生了一种线性变换，这个矩阵就代表了这种线性变换。
如第一中Loop Interchange是交换两层循环的位置，相当于循环迭代空间旋转90度，所以对应的变换矩阵是$[(0,1),(1,0)]$,
后面的变换也都类似，通过画图就可以比较清楚看出来所用的矩阵应该是什么。
么模变换的理论就是如果我们采用了一个线性的循环变换，对于的矩阵是A,那么对于每个数据依赖关系v,变换后数据依赖关系就变成Av了。
我们知道，对于一个循环，初始所有的数据依赖关系都是非负的。同样，对于一个合法的变换矩阵A,它必须使得对于每个数据依赖关系v,
Av也必须是一个非负的向量。这个就是判断一个循环变换是否合法的基本准则。

mathe

上面已经给出了如何判断一个变化是否合法的准则，现在来看看这些变换对应的矩阵有些什么特点。
首先，显然这些矩阵都必须是整数构成的，这样变换以后才能够保持循环变量都是整数。
其次，注意一下可以看出，每个基本变换的矩阵对应的行列式值都是1或-1.我们知道矩阵乘积的行列式等于行列式的乘积，
而两个变换复合对应的矩阵正好是两个变换对应矩阵的乘积，我们可以知道，所以这些变换以及它们的复合变换所对应的矩阵的
行列式值都是1。这样的矩阵被称为么模矩阵，这也是为什么这些循环变换会被称为么模变换。
有了这个理论基础以后，我们可以得出，像上面这样的一种3个变换的复合过程中，很有可能单单采用3个变换中任何一个都不是不合法的，
但是它们的复合变换确是合法的，我们需要做的只需要判断在这个复合变换对应的么模矩阵的作用下，是不是所有的依赖向量还是非负的就可以了。

下面来讨论如何将么模变换应用于并行化，很奇怪我看过的资料都没有讨论这方面的内容。也有可能么模变换功能还是太弱了些，所以并行化
机会还并不是很多，所以很少有人讨论。最近Stanford有些人研究一种么模变换的推广，叫做仿射变化，么模变换是对每个循环中所有语句采用了一种相同的
线性变换（么模变换），而仿射变换将对循环中每个语句都可以采用一个不同的线性变换，从而有更大的变换自由度，比如还可以表示Statement Shift,
Loop Fusion等变换，所有有更多的并行化机会。
只是仿射变换理论还不很成熟，至少我看了Amy W. Lim的毕业论文（讨论如何用仿射变换进行并行化)后发现其中有两个比较大的错误。但是如果手工采用这个变换
应该可以经常发现很多很好的机会。关于这方面资料，可以在Stanford的suif编译器的网页上找到(http://suif.stanford.edu/papers/ ,查找Affine Partition),而且SUIF编译器（研究性的编译器）已经将仿射变换加以实现。
而Amy W. Lim的毕业论文讨论的更加详细，在他在stanford的个人主页上应该可以找到（不过我刚才通过google没有搜索到）。
现在让我们回归主题，继续讨论如何将么模变换用于并行化。我们假设每个依赖向量都可以写成一个K维列向量 （循环总共K层），总共有S个依赖向量。
我们可以将这S个依赖向量放在一起形成一个K*S维矩阵D。如果现在有某个么模矩阵T使得TD (变换以后的依赖向量构成的矩阵）某一行全部是0，
那么也就是说，经过变换以后，这一行对应的循环中没有跨循环的依赖关系，也就是说，变换以后这一层循环对应的所有迭代都可以并行执行。
如果这一层是循环的最里层，那么变换以后的代码是一个最里层所有代码可以并行的程序，通常，这种程序我们可以:
i)将循环展开，让所有指令在指令级自动并行
ii)通过使用SIMD指令(Single Instruction Multiple Data)，比如SSE指令，让这个循环不同迭代中的代码并行执行
如果这一层是最外层循环，由于并行颗粒度通常比较大，我们可以采用多线程来并行。

如果阅读一下Amy W. Lim的毕业论文，我们可以发现，就是将变换扩展到仿射变换的范围，我们还是可以非常容易的找到让代码极大并行化的线性变换（只是仿射变换需要的循环变化方案非常不成熟，论文上提到的方法有问题）。

无心人

看着这么多东西
一个字
累

俺是菜鸟，俺不懂
还不如自己直接写汇编呢
（水平有限，嘿嘿）

gxqcn

楼上的，
这好比是一个讲堂，你听不懂或不爱听，尽可以打瞌睡，或悄悄离开，
犯不着这么乱嚷嚷，影响他人学习。。。

mathe

编译器本来就不是个简单的东西，看不懂也很正常。
有些东西，不是简单用汇编就可以写出来的。我们都知道，如果使用的算法不对，一个程序，你就是再怎么用汇编优化都没有用。

lyg_wangyushi

学习了，越来越崇拜mathe了，希望以后还能有这样的讲解，虽然不是很懂，不过受益非浅

mathe

现在举几个关于仿射变换的例子，比如代码：

1. for(i=0;i<N;i++){
   
2.     for(j=0;j<N;j++){
   
3.        A[i]+=...;
   
4.        ... = A[i+1];
   
5.     }
   
6. }

复制代码

对于这样的一个代码，它有一个关于循环I的跨循环的逆依赖关系。可以看出，
对于这个代码，它的数据依赖关系也是线性的，但是更复杂，不能够用依赖向量来表示.
我们现在看看仿射变换怎么处理这段代码。
通常仿射变换只处理循环中若干层，比如这里我们可以只重构最外层的循环。
那么我们重建一个新的循环变量I, 对于每个语句，它都是旧的循环变量的线性组合。
比如对于第一个语句，我们取
  I=1*i+0*j=i
对于第二个语句，我们取
  I=1*i+0*j+1=i+1
对应代码经过变换后成为:

1. I=0;
   
2. for(j=0;j<n;j++)
   
3.    A[I]+=...
   
4. for(I=1;I<n;i++){
   
5.    for(j=0;j<n;j++){
   
6.        ...=A[I];
   
7.        A[I] +=...;  
   
8.    }
   
9. }
   
10. I=n-1;
    
11. for(j=0;j<n;j++)
    
12.    ...=A[I];

复制代码

可以看出来，这两份代码中，第一个代码关于最外层循环是有跨循环的数据依赖关系的，
所以我们无法对最外层循环做并行化;但是第二个程序不同，它关于最外层的循环的所有循环迭代
之间没有任何数据依赖关系，如果循环体足够复杂，那是非常适合转化成多线程，让多个核同时运行，
从而可以提高在多核上的速度。上面这种变换编译器中称为statement shifting（就是平移了一个语句）

下面举一个更加复杂一点的例子

1. for(i=i0;i<n;i++)
   
2. for(j=j0;j<n;j++){
   
3.      A[i][j]+=B[i-1][j];                    //(s1)
   
4.      B[i][j]=A[i][j-1]*B[i][j];             //(s2)
   
5. }

复制代码

可以如图表示出所有的数据依赖关系

对于这个代码，如果使用么模变换，那么所有同一个循环迭代内部的(s1),(s2)由于被统一处理，
我们可以看成它们两被合并后变成一个点，这样我们得到的将是一副连通的图，也就是数据依赖
关系将会让所有的计算依次依赖，所以我们必须串行执行代码
但是换成仿射变换，由于每个语句可以有不同的线性变换，所以同一循环迭代内部的(s1),(s2)这时候
是两个不同的点，如上图，我们可以将整个被执行的数据划分成一些相互不连通的子图。也就是不连通的
部分实际上相互之间没有数据依赖关系，所以我们可以放在不同的核上并行处理
对于这个变换，实际上的操作是引入一个新的外层循环变量P,让
P针对S1是i-j,P针对S2是i-j+1
为了方便的实现这个变换，我们可以首先直接在外层添加一个多余的P循环，代码变成

1. for(P = 1-n;P<=n;P++){     
   
2.   for(i = 1;i<=n;i++){
   
3.      for(j = 1;j<=n;j++){
   
4.         if(P == i-j)
   
5.            A [i][i- P] += B[i-1][i-P];
   
6.         if( P == i-j+1)
   
7.            B [ i][i-P+1] = A[i][i-P]*B[i][i-P+1];
   
8.      }
   
9.   }
   
10. }

复制代码

下一步，我们可以从每个语句附带的限制条件解出j仅仅在某些特殊情况才会被执行，然后优化掉内层循环就可以（很多编译器本身就支持这种优化）
得到：

1. for(P = 1-n;P<=n;P++){
   
2. if(P>=1)
   
3.          B[P][1]=A[P][0]*B[P][1];
   
4.   for(i=MAX(1,P+1) ;i<=MIN(n,P+n-1);i++){
   
5.         A [i][i- P] = A[i][i-P]+B[i-1][i-P];
   
6.         B [i][i-P+1] = A[i][i-P]*B[i][i-P+1];
   
7.    }
   
8.    if(P<=0)
   
9.         A[i][i-P]=A[i][i-P]+B[i-1][i-P];
   
10. }

复制代码

变化以后的循环的最外层就可以完全并行化了。
上面的变化虽然导致最外层可以并行化，但是如果实际运行，我们会发现性能反而会变差，这个是因为现在的计算机的内存访问模式都使用了多级缓存(cache)
结构，对于原始的代码，对于内存的访问都是顺序进行的，而变换以后的代码，内存访问顺序就乱序了，这个会导致性能极大的降低，
那么我们有没有办法兼顾呢？编译器里面有另外一种技术，称为循环分块（Blocking，不过可以这种循环变换无法集成到仿射变换中），
我们可以先对并行化后的代码进行循环分块（单独对P循环分块就可以了），然后对于每个分块，重新反变换回来就可以了（或者采用另外一种仿射变换
使得内存访问在变换后最优）
变换的效果相当于对于原循环迭代空间入下图进行划分，不同颜色的迭代划分到不同线程（当然假设n很大）

其中每个线程划分到了一段P.
变换后每个线程的代码大概类似,假设对于对应的线程，划分到的P的范围是从MINP到MAXP

1. for(i = 1;i<=n;i++){
   
2.    if(i-MAXP>=1)
   
3.            A[i][i-MAXP] = A[i][i-MAXP]+B[i-1][i-MAXP];
   
4.    for(j = MAX(1,i-MAXP);j<=MIN(n,i-MINP)){
   
5.            A[i][j] = A[i][j]+B[i-1][j];
   
6.            B[i][j] = A[i][j-1]*B[i][j];
   
7.    }
   
8.    if(i-MINP+1<=n)
   
9.            B[i][i-MINP+1]=A[i][i-MINP]*B[i][i-MINP+1];
   
10. }

复制代码

mathe

关于循环优化的理论比较复杂，我们现在暂时切换到一些简单些的优化吧
我们经常在写代码时会使用递归调用，我们知道递归调用通常代价比较大，
但是递归调用比较符合我们的思考方式，写代码比较容易。
其实实际上对于一种比较特殊的递归调用，编译器可以直接帮助我们做优化，
将递归调用直接优化为非递归调用，这一类递归调用在编译器里面称为尾递归调用
(Tail Recursive Call).
比如函数

1. int primes(int a, int b){
   
2.     if(a==0)return b==1;
   
3.     if(b==0)return a==1;
   
4.     return primes(b,a%b);
   
5. }

复制代码

像上面的一个代码，只有在代码的最后一行调用函数自身的代码称为尾递归，大部分现代
的编译器都能够对上面的代码进行优化，将递归优化成一个循环，因为其做法只需要修改所有
函数输入值，然后跳到函数入口就可以了。
有些尾递归调用可能不是很容易看出来，比如

1. int primes(int a, int b){
   
2.     if(a==0)return b==1;
   
3.     if(b==0)return a==1;
   
4.     if(a<b){
   
5.         return primes(a,b%a);
   
6.     }else{
   
7.         return primes(a%b,b);
   
8.     }
   
9. }

复制代码

虽然上面代码有两次递归调用，但是两次都是尾递归，它们分别在两个不同分支的最后面。
此外另外一种比较常见的尾递归调用形式是:

1.    double honic(int n){
   
2.         if(n==1)return 1.0;
   
3.         else return 1.0/n+honic(n-1);
   
4.    }

复制代码

这里，递归调用在一个表达式中间，虽然它是最后一次计算(+操作)的一个参数，但是由于加法
满足结合率，我们可以改变所有累加的顺序，所以上面的编译器也可以自动转化为循环。
同样，计算

1.    double honic(int n){
   
2.         if(n==1)return 1.0;
   
3.         else return honic(n-1)+1.0/(n*n);
   
4.    }

复制代码

编译器也可以自动优化，（虽然后面还有计算1.0/(n*n),但是编译器可以自动优化，加法满足交换率).
像下面代码

1.    double Fib(int n){
   
2.         if(n<=1)return 1.0;
   
3.         else return Fib(n-1)+Fib(n-2);
   
4.    }

复制代码

这个代码中，两次递归调用中，只有后面那次调用(Fib(n-2))可以转化为循环的一部分，但是前面那次递归调用无法转化。
所以对这样的代码，最好还是通过手工优化。（其实，这里的函数Fib不含有任何指针之类会引起歧异的代码，从理论上，
编译也应该能够优化掉，但是实际上现在通常的编译器不会做这样的优化）

gxqcn

曾在CSDN上看到有人提关于“尾递归”的问题，当时只能从字面上理解，
听了楼主的讲解，才有了更深刻的认识，谢谢！