Chinaunix首页 | 论坛 | 博客
  • 博客访问: 1994021
  • 博文数量: 960
  • 博客积分: 52560
  • 博客等级: 大将
  • 技术积分: 13131
  • 用 户 组: 普通用户
  • 注册时间: 2008-07-31 14:15
文章分类

全部博文(960)

文章存档

2011年(1)

2008年(959)

我的朋友

分类: C/C++

2008-08-01 17:11:42

下载本文示例代码
简介

在《C Templates: The Complete Guide》一书中(以下简称书),提出了模板元编程最早的实际应用之一:在数值运算中进行解循环优化。
而本文的标题是噱头!本文的真正目的是指出这种优化措施在增加复杂性的同时,并不一定能明显改善效率。应当谨慎使用该技术——默认不使用该技术,在点积计算确实是效率瓶颈时考虑采用该技术,并认真测试该技术是否真能提高效率。

背景

数值运算库中常常需要提供向量点积(dot_product)运算。其定义用C 代码描述也许更清楚~

template

T dot_product(int dim,const T v1[],const T v2[]) {

	T result=0;

	for (int i=0;i我们可以使用这个函数,求2个向量的点积
int v1[] = {1,2,3};

int v2[] = {4,5,6};

int r1 = dot_product(3,v1,v2);

得到r1=32

书中指出:“这个结果是正确的,但是在要求高效率的应用中,它耗费了太多的时间”,对于这类特殊的问题,“简单的把循环展开”,如:r1=v1[0]*v2[0] v1[1]*v2[1] v1[2]*v2[2]
“反而会好得多”。

如何便捷的展开循环?将每个dot_product(dim,...) 手工重写成展开式?还是设计一整套函数: dot_product_dim_2,dot_product_dim_3,... ?
无疑这是一个冗长乏味易出错的方案。

书中提出了一种使用模板元编程解决该问题的方案,让我们来看看他是怎么做的。

模板元编程——解循环

// 首先是递归模板

template

struct DotProduct {

	static T execute(const T v1[],const T v2[]);

}

template

T DotProduct::execute(const T v1[],const T v2[]) {

	return v1[0]*v2[0]   DotProduct::execute(v1 1,v2 1);

}



// 递归边界模板

template

struct DotProduct<1,T> {

	static T execute(const T v1[],const T v2[]);

}



template

T DotProduct<1,T>::execute(const T v1[],const T v2[]) {

	return v1[0]*v2[0];

}

我们可以使用DotProduct 来进行点积计算了:
int v1[] = {1,2,3}; int v2[] = {4,5,6};

int r2 = DotProduct<3,int>::execute(v1,v2);


计算r2的函数,在DotProduct<3,int>::execute 被实例化时就确定了。
编译器将会实例化
int DotProduct<3,int>::execute(const int v1[],const int v2[]) {

	return v1[0]*v2[0]   DotProduct<2,int>::execute(v1 1,v2 1);

}

然后编译器继续实例化
int DotProduct<2,int>::execute(const int v1[],const int v2[]) {

	return v1[0]*v2[0]   DotProduct<1,int>::execute(v1 1,v2 1);

}

这里,我们有一个 DotProduct<1,T> 的偏特化版本,计算函数将被实例化为
int DotProduct<1,int>::execute(const int v1[],const int v2[]) {

	return v1[0]*v2[0];

}



而这3个函数都足够短小,编译器通常会为它们做inline工作,使得r2的计算被展开为
r2 = v1[0]*v2[0] v1[1]*v2[1] v1[2]*v2[2];

DotProduct<3,int>::execute(v1,v2) 的调用语法不是很友好。
书中还设计了一个helper函数。
template

T dot_product(const T v1[],const T v2[]) {

	return DotProduct::execute(v1,v2);

}

如同STL的make_pair、bind2nd,这个模板函数将模板参数T的推导工作交给编译器,使DotProduct的界面更友好,更不容易出错。当然,维度DIM必须指定。

现在可以按如下方式使用 :

int r2 = dot_product<3>(v1,v2); // 是不是和 int r1 = dot_product(3,v1,v2); 十分相似?

本文以下部分将从编译器生成的机器码层次来验证这一技术,同时与普通循环点积计算做比较。
测试环境 :Windows Xp sp2,Microsoft Visual Studio 2005 Team Suite,release默认设置

验证1:

该技术真的能和预想的一样,实例化多个函数,并将它们inline,达到解循环的效果吗?
让我们写个短小程序验证一下:

int main() {

	int v1[] = {1,2,3};

	int v2[] = {4,5,6};

	int r1 = dot_product(3,v1,v2);		//循环版本

	int r2 = dot_product<3>(v1,v2);	//解循环版本

	cout<


程序输出当然是2个32,让我们进入反汇编看看实际工作情况:
令人吃惊的是,前4行代码根本没有生成任何机器码
生成机器码的语句只有2条输出语句和return 0;
cout< >::operator<< (40203Ch)]

mov	ecx,eax

call	dword ptr [__imp_std::basic_ostream >::operator<< (402038h)]

这些指令的的含义是什么呢?

我们知道 cout< (cout< cout.opeator<<(...) 返回的是自身引用,所以可以进行链式输出

前4条指令执行后:
堆栈应该是:
std::endl
20h
ECX 保存的是 std::cout

std::cout 的类型是std::basic_ostream >
所以第5条指令正是对cout的成员函数“ opeator<<”进行调用
(this指针 std::cout,已经存入ECX,而且可以猜想这个成员函数正是 operator(int) 重载)
而参数,对的,是call指令前,堆栈栈顶“20h”—— 32的16进制。
第5条指令将执行 cout.operator<<(20h); !

让我们把 r1被编译器偷偷替换成20h的事放在一边,继续解析后2条指令:
cout.operator<< 使用__thiscall调用约定,函数自己弹出堆栈中参数20h, 返回自身引用,返回值存入EAX。
所以第5条指令执行后:
堆栈应该是:
std::endl;
EAX 保存的是返回值,即是std::cout

第6条指令 mov ecx,eax 将返回值重新放入ECX,设置好this指针。
第7条指令 继续调用 operator<< 的另一个重载版本。(这时 std::endl 在栈顶)
完成语句 cout.operator<<(endl);

cout<

这个原本就短小的程序,被编译器“篡改”为更短小的程序:

int main() {

	cout<<32<如果改用 printf 输出:
printf(“%d\n”,r1); printf(“%d\n”,r2);
从反汇编可以更清楚的看出,编译器实际做出的代码是:
printf(“%d\n”,0x20); printf(“%d\n”,0x20);

(详细代码见 sample1)

比较合理的解释是,编译器知道了太多的上下文,做出了足够的优化,将本来应该运行时得出的结果—— 一个不变的结果 ——放入的可执行文件中,直接输出。

为了验证模板元解循环技术,我们要另外想办法。

验证2:

我们需要“阻挠”编译器,让它不能猜测出运行结果,从而生成真正的计算代码。
这里使用了一个比较简单的方案:

template

void __stdcall MyGenerate(OutIt first,OutIt last,const char *name) {

	cout<<"generating "< oit(cout," ");

	copy(first,last,oit);

	cout<

  

用clock 函数返回值填充一个序列,即使编译器“胆大包天”,也不敢猜测clock返回结果了~
多余的 name 参数和一些输出语句是因为在release模式下,watch窗口中看不到值。
而__stdcall 是让 MyGenerate后不生成多余的平衡堆栈代码。

按照如下方式使用 :

MyGenerate(v1,v1 3,”v1”);

MyGenerate(v2,v2 3,”v2”);

int r1 = dot_product(3,v1,v2);

MyGenerate(v1,v1 3,”v1”);

MyGenerate(v2,v2 3,”v2”);

int r2 = dot_product<3>(v1,v2);

cout<

  

仍然不够, 编译器会将r2的计算,推迟到 cout<

所以我们再加入一个函数

void ForceCalcResult(int result) { }

强制编译器立刻计算点积。
同时,使用函数指针调用,避免编译器对该函数inline。

程序主干如下 :

int main() {

	const int dim = 3;

	int v1[dim];

	int v2[dim];

void ( *const ForceCalcResult_non_inline)(int) = ForceCalcResult;



	MyGenerate(v1,v1 dim,"v1");

	MyGenerate(v2,v2 dim,"v2");

	int r1 = dot_product(dim,v1,v2);

	ForceCalcResult_non_inline(r1);



	MyGenerate(v1,v1 dim,"v1");

	MyGenerate(v2,v2 dim,"v2");

	int r2 = dot_product(v1,v2);

	ForceCalcResult_non_inline(r2);



	cout<这样,编译器就屈服了,乖乖的把r1和r2的计算代码展现出来。
//MyGenerate(v1,v1 dim,"v1");

mov		eax,offset string "v1" (402134h) 

lea		edi,[esp 24h] 

lea		ebx,[esp 18h] 

call		MyGenerate (401220h) 

//MyGenerate(v2,v2 dim,"v2");

mov		eax,offset string "v2" (402138h) 

lea		edi,[esp 18h] 

lea		ebx,[esp 0Ch] 

call		MyGenerate (401220h)

从这4条指令,我们可以看出v1和v2的地址:
v1[0]:esp 18h, v1[1]:esp 1Ch, v1[2]:esp 20h, v1[dim]:esp 24h

v2[0]:esp 0Ch, v2[1]:esp 10h, v2[2]:esp 14h, v2[dim]:esp 18h

让我们先看模板解循环版本:
//int r2 = dot_product(v1,v2);

mov		edi,dword ptr [esp 10h] 

mov		edx,dword ptr [esp 14h] 

imul		edi,dword ptr [esp 1Ch] 

imul		edx,dword ptr [esp 20h] 

// edi = v2[1]*v1[1]

// edx = v2[2]*v1[2]



mov		eax,dword ptr [esp 0Ch] 

imul		eax,dword ptr [esp 18h] 

// eax = v2[0]*v1[0]

add		edi,edx 

add		edi,eax 

// edi = edi edx eax = v2[1]*v1[1] v2[2]*v1[2] v2[0]*v1[0]

循环被解开了!利用3个寄存器edi,edx,eax,最终结果应该是保存在edi中。
再看接下来的代码

//ForceCalcResult_non_inline(r2);

push		edi

call		ForceCalcResult (401000h)

结果确实是存放在edi中的。
我们接着看普通“循环”版本:
//int r1 = dot_product(dim,v1,v2);

mov		esi,dword ptr [esp 10h] 

mov		eax,dword ptr [esp 14h] 

imul		esi,dword ptr [esp 1Ch] 

imul		eax,dword ptr [esp 20h] 

// esi = v2[1]*v1[1]

// eax = v2[2]*v1[2]

mov		ecx,dword ptr [esp 0Ch] 

imul		ecx,dword ptr [esp 18h] 

// ecx = v2[0]*v1[0]

add		esi,eax 

add		esi,ecx 

// esi = esi eax ecx = v2[1]*v1[1]   v2[2]*v1[2]   v2[0]*v1[0]



//ForceCalcResult_non_inline(r1);

push		esi  

call		ForceCalcResult (401000h)

几乎相同的代码,使用的是esi,ecx,eax,结果保存在esi中。
循环同样被解开了!

编译器在我们的重重算计下,被迫在运行时计算结果,同时将运算方法展现出来;但同时,它依然不屈不饶的向我们展示它的强大优化能力!

(详细代码见sample2)

验证2 :

在验证2中,编译器知道了 const int dim = 3; 这一上下文。从而将普通循环版本的点积函数进行展开。

让我们来看看,dim取更大值时,编译器会如何。

dim=10:
普通“循环”点积计算被展开,使用esi,eax,ecx,edx,10条mov与imul指令,9条add指令,最终将计算结果存入esi
模板元点击计算也被展开,使用edi,eax,ecx,edx,10条mov与imul指令,9条add指令,最终将计算结果存入edi

dim=11:
循环点积计算发生了变化,代码如下:

//MyGenerate(v1,v1 dim,"v1");

00401017  mov	eax,offset string "v1" (402134h) 

0040101C  lea	edi,[esp 68h] 

00401020  lea	ebx,[esp 3Ch] 

00401024  call	MyGenerate (401280h) 



//MyGenerate(v2,v2 dim,"v2");

00401029  mov	eax,offset string "v2" (402138h) 

0040102E  lea	edi,[esp 3Ch] 

00401032  lea	ebx,[esp 10h] 

00401036  call	MyGenerate (401280h) 



// v1[0]:esp 3Ch, v2[0]:esp 10h

//int r1 = dot_product(dim,v1,v2);



0040103B  xor	ebp,ebp 

// int result = ebp=0;

0040103D  xor	eax,eax 

// eax=0;

0040103F  nop

// align



// i=eax;

// loops:

00401040  mov	ecx,dword ptr [esp eax 10h]

// ecx = *(v2 i);

00401044  imul	ecx,dword ptr [esp eax 3Ch]

// ecx *= *(v1 i);



00401049  add	eax,4

//   i;



0040104C  add	ebp,ecx 

// result =ecx



0040104E  cmp	eax,2Ch 

// if (i<11)



00401051  jl		main 30h (401040h) 

// goto loops;

编译器已经不能“忍受”这样长度的循环展开,将其老老实实的翻译成真正的循环。
但是,编译器仍然对函数做了inline工作。

对于模板元解循环,因为代码要求11层函数调用,编译器能做的只有将11层调用inline,得到的是就地展开的计算。

dim=12:

循环点积计算: 编译器甚至连inline也不做了,进行正式的函数调用。
模板元点击计算:依然就地展开。

dim=32:

循环点积计算:不展开,调用
模板元点积计算:
比较有趣的是,编译器也没有将函数就地展开,而是分成3步。就地计算前9个乘法然后与DotProduct<23,int>的返回值相加。DP<23>也没有就地计算全部,而是计算前11个,然后与DP<12>的返回值相加。
3步计算当中都没有循环。

值得注意的是,循环点积计算函数也没有很老实的按照源代码的方式进行计算,而是进行了不完全的循环展开,类似于如下代码:

template< typename T>

T dot_product(int DIM,const T v1[],const T v2[]) {

	const int step = complierKnow;

	T results[step] = {0};

	int i=0;

	for (;i step-1

DIM和step似乎没有什么固定的规律。

(详细代码见sample2,取不同的dim值即可。)

验证3:

sample2中,编译器对被计算的向量的上下文依然知情:v1,v2在main函数的栈中。

sample3做了小小的改动:dot_product_loop和dot_product_unloop将无从得知传递给它们的向量在何处。(当然,在main函数中,仍然使用函数指针来调用这2个函数,使得编译器不做inline)

sample3得到的结果和sample2基本相同,只是对向量寻址由esp v1_offset,esp v2_offset变为[esp 4],[esp 8]

(详细代码见sampl3)

验证4:

在sample3的基础上,通过控制台输入读取dim,使得编译对其不知情。
当然,在这种情况下,是无法使用模板元解循环的,因为它要求dim是编译时常量。

在这种情况下,循环点积计算终于被编译器翻译成真正的“循环”了。

总结:

循环版本有更大的灵活性:在dim为编译时常量时,编译器会根据其大小,进行完全解循环或部分解循环。同时也支持dim为运行时的值。
模板元版本必须使用编译时常量作为dim,而且总是完全解开循环。

循环版本与模板元版本最终都会使用相同次数的乘法与运算,区别在于乘法指令和跳转指令的数目。
模板元版本在dim较大的时候,可能会使用跳转指令,将部分工作交给更小维度的点积计算。但是乘法指令数目与维数一样。
循环版本可能会使用更多的跳转指令,使得乘法指令数目大大减少。

多出的跳转指令会占用运行时间,但也能减少目标代码体积。权衡应该使用那一种版本与权衡某个函数是否应该inline十分相似。
书中17章最后部分也提醒读者,不计后果的解循环并不总是能优化运行时间。

借用大师的观点来结束本文 —— 优化的第一原则是:不要优化。优化的第二原则(仅适用于专家)是:还是不要优化。再三测试,而后优化。
《C 编程规范》 第8条

下载本文示例代码
阅读(606) | 评论(0) | 转发(0) |
给主人留下些什么吧!~~