如何计算一张32位图片的储存?复数相乘的计算后为什么要除以8?…

自古计算圆周率是数学界一项流行运动,各大数学家争相破记录以名垂青史。想象有人为圆周率,算的不是圆周率而是寂寞啊!

自有圆周率,计算比的是数学;后有现代数学,计算比的是寂寞;自从有了计算机,计算变成程序员们(另一种)练手的健康活动:锻炼编程技术之余可比肩历史伟人,看官也来一发吧!不懂数学没关系,计算圆周率并不需要多少数学知识,你需要的只是一些编程基础、一点数据结构知识,还有本文准备给你的:算法。

(题外话:如果本文中所有公式都显示成叉,说明你在使用古董级浏览器,请升级……)

1)本文介绍编程计算圆周率的主流算法,从简单到复杂任君选择,目的是让读者看完能动手算一把。
2)不涉及具体的编程语言技术,也不涉及数学证明。
3)所有算法都不限编程语言,你可以用你最擅长的语言来实现。不过在科学计算领域,C/C++兼顾性能和便捷,依然是首选,文中少量示例都用C/C++。
4)即使你没时间动手,本文也可带你浏览一下科学计算王国。
【读懂本文,你需要一点点的……】
1)编程基础
2)数据结构基础
3)高等数学基础

千言万语不如先来一发:

请把以上代码拷进C语言开发环境里运行,瞬间出结果如下:
pi=3.79

恭喜!14位小数无一错误!你已经破了公元1400年的纪录!祖冲之那7位小数真的弱爆啦!……只是开了个玩笑,算圆周率前请先拜大神,可保佑不溢出:
祖冲之( 公元429年4月20日─公元500年),我国杰出的数学家,科学家。曾计算出圆周率7位小数,是之后近1000年的世界纪录保持者。

凡是讨论算法,必谈时间复杂度。祖冲之的年代,圆周率源自几何也算以几何,原理是在圆周割成多边形来计算周长,称为几何算法,时间复杂度高、计算量极大。经过现代数学和计算科学的发展,圆周率计算方法变得非常高效,例如上面的C语言例子用的是以下无穷三次级数:
不用计算机,笔算也能算出好几位,时间复杂度是O(10^(N/3*2))(N是十进制位数,下同),但仍不足以计算成千上万位。另外,例子中用double(双精度浮点数)类型来计算圆周率,但编程语言支持的浮点类型最多就十几、二十位几小数,显然上面那种简单的程序无法算出更精确的pi。为了让计算机处理超长的小数,我们要用很大的数组来储存小数(称为【高精度计算】),并为此实现专门的加减乘除、开方算法。本文后续将详细介绍高性能的圆周率计算方法和与之配套的高性能高精度计算。

圆周率计算方法极多,但不是每种都能高效实用。现代圆周率计算方法主要有两种:

  1. 基于无穷幂级数(简单,O(N^2))
  2. 基于超长小数运算(高速,O(N logN))


第一种方法基于pi可以***成一组幂级数和,例如等一下我会介绍的贝利-波尔温-普劳夫公式:
非常简单有效的算法,用最基本的高精度四则运算即可实现,没有高精度的分母,没有两个大数相乘,更没有开方,时间复杂度是O(N^2) ,推荐从这个上手。

第二种方法包括稍后介绍的高斯-勒让德算法,依靠上世纪60年代发明的神器快速傅里叶变换而达到时间复杂度O(N(logN)^2),位数多时优势明显,是世界纪录和众多圆周率计算软件的算法,包括著名的SuperPI,但用到的高精度计算技术稍复杂。

在圆周率计算和各种科学计算中,都用数组来储存一个很长的小数(或者超大的整数,或小数点前后都有很多位)。例如,下面表示用【高精度数组】a储存的圆周率,数组每个单元保存8位十进制小数:
等同于:

其中是进制。一般地,小数A用高精度数组a基于d进制储存可表示为:
计算机的储存空间有限,所以n是一个有限的整数。进制d一般是10或者16的整数次幂,其中16整数次幂的方法对储存的利用率更高、进位计算更快速,但输出时转成十进制需要很大的计算量(O(N^2)),有时甚至比计算pi本身更耗时间。

对高精度数组表示的小数进行各种加减乘除运算,称为【高精度计算】。例如,最简单的高精度计算是加法,把两个高精度数组对位相加,即可得出两个小数的和高精度小数:
代码看起来是这个样子的:

如果结果中若干项超过进制d,则对其进行进位运算。C语言进位的代码看起来是这个样子的:

进位可能也造成a[k-1]也超过了d,于是进位从最低位a[n]一直到最高位a[0]反向遍历执行。

减法(或者负数加法)的原理也类似,对位相减后执行退位操作(略)。

【高精度乘法】基本计算方法是【卷积】,小学时教我们的竖式乘法本质也是卷积。所谓卷积,就是多项式每一项都与另一个多项式的每一项相乘、累加:
小学教的竖式乘法就是卷积,被乘数和乘数各位两两相乘、相加、进位:

高精度数组a和b相乘的数学表达是:
所以C的每一项是:

上面的式子看傻了没关系,高等数学就是喜欢把简单的东西写得很复杂……其实原理极简单,代码看起来是这个样子的:

同样要进行进位操作。两个n位精度的小数相乘,虽然结果有2n-1项,但其实精度只有前面n项是准确的,后n-1项丢弃或不予计算。卷积的时间复杂度是O(nm)。

如果加减乘对你来说都没有压力,那【高精度除法】也许能给你带来一些麻烦。高精度除法基本算法就是小学教的竖式除法:试除一位x、计算【被除数-(x×除数)】得余数,余数进行退位运算后作为被除数重复上面的过程。如果除数是一个整数,这个并不会很困难,试除用的是被除数的前几位除以除数。但如果除数也是高精度数组就有点麻烦了,除数也要用前几位作为试除位来试除。这种基本除法算法时间复杂度是O(nm),n是结果精度,m是除数精度,与被除数精度无关。实际上这种代码实现起来比后面的高性能计算方法更长更麻烦(因为后面介绍的方法把除法转化成乘法,可直接调用乘法的代码),所以这里略去。

在你实现高精度运算时,需要决定数制,或写出兼容任何数制的代码(但这并不好)。所谓数制,就是我们平时使用的十进制、十六进制。唯一不同的,高精度的数制不是其显示出来的数制,而是其每个单元储存数据的上限。例如规定高精度10000进制,表示一个单元达到或超过10000时,应该减去10000,并对其高一位的单元加1(进位),或是一个单元小于0,则加上10000,其高一位减1(退位)。

一般我们选用10或16的整数次幂作为高精度的数制,两者各有优缺点。10整数次幂的数制运算完后无需转数制可直接输出十进制的结果;而16整数次幂进制进位速度高,只需要用编程语言的位运算符,性能远好于取模、除法运算,而且下面的BBP公式和十六进制是绝配,通用所有数制的程序几乎不可能利用有这些优点,我们没必要写出通用的程序。

因此,有时我们可能选择用16整数次幂进制计算,然后转成十进制输出。转换算法很简单,对于任意一个高精度16整数次幂进制的小数,重复下面的步骤:

  1. 整数部分直接用十进制输出
  2. 整个高精度数组乘以10^n,并用原来的数制执行进位
  3. 回到1重复这个过程,直到输出十进制位数达到十六进制位数的1.2倍()


以上计算每一个循环能输出n位,总时间复杂度是O(N^2):这个是个大问题,超过了高斯-勒让德算法的时间复杂度,所以如果用16整数次幂进制跑高斯-勒让德算法然后再转十进制会得不偿失,转换的速度甚至比计算更慢。幸好,这个计算易于实现多线程并行运算,很容易充分发挥CPU超过4个线程的性能。

现在你已经学会了基本的高精度计算,用BBP来算已经难不倒你了。

BBP公式是圆周率计算的众多无穷幂级数方法之一,均具有O(N^2)的时间复杂度,优点是实现简单、极容易实现多线程并行运算、内存使用效率高,使用最为广泛。


其中每个方块都表示一个高精度数。留意左边1/16^k越来越少,总有某一项之后小于你设置的高精度数组的精度,之后就不用计算了,可输出结果。计算过程只有高精度加法、减法,高精度数乘除非高精度整数。计算总项数n正比于计算精度,每轮计算涉及的高精度计算时间复杂度都是O(N),所以总复杂度是O(N^2) 。

刚才说,BBP和十六进制是绝配,不只是因为免去左边除以16的O(N)除法,还因为BBP公式在十六进制计算时可以只计算其中某些位而无需计算前面的未,例如从n位计算到n+m位,时间复杂度是O(n*m)。基于这个特点,BBP不是圆周率幂函数计算公式最高效的一个,却最适合多线程并行和分布式计算,还能用来验算其它算法算出来的十六进制结果。考虑到一点跑题,不在这里详述,有兴趣可看看(较简单)。

BBP之类的O(N^2)方法很简单,但在世界纪录级别(万亿位)难等大雅之堂,我们需要的一些O(N logN)级别的神器。其中最简单的是GLA,SuperPi用的就是这个。

GLA写出来很简单。
初始化:
迭代:
最后计算pi为:

其中,除了p外,其余(a,b,t,pi)都为高精度数,例如计算100万位圆周率,则这4个高精度数都是100万位。迭代过程中,a和b会不断接近,迭代结束条件是a和b的差距超过其本身精度的一半(例如计算100万位,那结果a和b第一个不同的位在50万位后即可结束迭代)。p是普通的整型。

每迭代一次,a和b的差的数量级增加一倍,所以迭代次数为O(log N),计算100万位的圆周率约需要迭代20次。

GLA难点和主要计算量在高精度乘法和开方,乘法和开方的性能决定GLA的性能。要使GLA的时间复杂度达到O(N (log N)^2),你需要下面的算法。

基本的高精度乘法用卷积,时间复杂度O(nm),用它的话我们GLA的O(N (log N)^2)就报销了。我们要祭出神器中的神器:FFT——快速傅里叶变换。有了FFT,对大数乘法、除法、开方通通有了过往不可想象超高性能,Pi的位数疯狂增长。

尽管我尽了最大努力,但本节仍是这篇文章最难读懂的部分,你可能会一头雾水——没办法,神器自有神器的架子,在我自行参透前网上几乎没找到清晰的讲解,所以写这部分结合了我自己的经验,着重解释最难理解的部分,并使用大量的图示,尽可能让过程更直观。尽管看起来很复杂,但代码核心部分也只有10-15行,因为其有着惊为天人的对称性,体现着让人叹为观止的数学美。

【傅里叶变换】和【快速傅里叶变换】:【傅里叶变换】本来不是用来算乘法的,而是用来***波(例如声音或者电磁波),而用来做卷积则是上世纪自有计算机后发现的。傅里叶变换将一组复数通过一系列计算转换成另外一组复数,转换前后复数个数一样,转换后的每个复数都和转换前每个复数相关。【快速傅里叶变换】泛指可以在O(N logN)时间复杂度内完成傅里叶变换的算法。根据卷积定理,两个多项式的卷积的傅里叶变换等于其各自的傅里叶变换对位相乘,如果结合快速傅立叶变换,就可以在O(N logN)的时间复杂度内完成高精度乘法。本文只介绍算法而不论数学原理和证明,其实我还没搞懂我会告诉你吗?

【算法总体流程】假设要对高精度数组A和高精度数组B以快速傅里叶变换计算卷积,总体流程如下:

  1. 把数组A转换成复数数组A';
  2. 对A'执行快速傅里叶变换,得复数数组A'';
  3. 同上方法,把高精度数组B转换成B',执行快速傅立叶变换得B'';
  4. A''和B''对应单元相乘,得复数数组C'';
  5. 对C''执行快速反傅里叶变换,得复数数组C';
  6. C'重新转换成高精度数组C,即为计算结果。

其中,转换就是普通的数据类型转换,整数转成一个虚部为0的复数,或者复数舍弃虚部,实部四舍五入成整数。
“快速反傅里叶变换”(IFFT)是“快速傅里叶变换”(FFT)的逆过程。事实上由于数学上的对称性,两者几乎是一样的,核心代码只用一套,无需分别写。

快速傅里叶变换有多种具体算法,以下介绍最简单的“库利-图基快速傅里叶变换”【C-T FFT】。

【转换】C-T FFT要求上图中所有复数数组统一大小,且为2的整数次幂,不小于相乘两个数的长度和,所以计算前需要先补0。假设计算1.11×2.22,高精度数组是十进制,1.11和2.22都是3个单元的精度,3+3=6,不小于6的2的整数次幂是8,所以整个流程看起来是这样的(用(a,b)表示复数a+bi,下同):

【转换精度问题】复数使用两个浮点数表示,计算机支持的64bit浮点数一般只有54bit的有效位数(GCC支持96bit long double类型并不通用)。假设高精度整数数组每个单元是 x bit,即不大于进制,FFT长度是,那么为了防止精度不足而造成计算错误,x和N必须满足:
这个条件相当苛刻,意味着 x 得不超过16,即不超过2^16=65536进制。如果用10整数次幂计算,不能超过10000进制。所以当你用的数制较大,转复数数组时必须拆开到65536进制或10000进制。当高精度数组尺寸达到约2^24时,还得进一步下降到256进制和100进制。

【结果精度问题】用n位d进制表示一个无穷小数时,误差不超过,那么两个小数相乘后,结果误差不超过(和加减法一样),仍然是n位,超过n位即使算出来也是不准确的,可以舍去。证明很简单,略。

【FFT 和 IFFT】刚才提到,傅里叶变换和反傅里叶变换有漂亮的对称性,IFFT只需在FFT前后再加一点东西:
其中,共轭就是复数数组中所有单元的虚部取负值(高中数学呀),归一化就是全部复数除以N(数组的大小)。由于我们计算高精度乘法时,输出结果虚部一定是0,所以最后一个共轭计算可以节省。

【FFT过程】公式神马的最讨厌了,要说清楚FFT,千言万语不如一张图:

即使你没看明白也一定能感受到这图的震撼,它有着各种完美的周期性和对称性。留意从顶上的输入数组到底下输出结果数组有很多连线,输入的任意单元和输出的任意单元之间有且只有一条通途。

其中W(k,N)是“旋转因子”,表示模为1,-k/N个圆周角的复数,即:

(*)所谓位元反转排序就是把数组单元的序数k的N位二进制全部位取反,例如:
4的二进制是100,把所有位反过来是001,是1的二进制,所以x4和x1交换。
N=8时,0、5等因为位元反转后仍是自身,所以x0、x2、x5、x7无需换位。

请细细观察图中的规律,直到你有点头绪了,看看下面我归纳的规律:
一、对于N=2^p的FFT,有p轮计算。
二、每一轮计算都两两组合。
三、决定如何组合两个单元的规律在于两个单元的序数的二进制,例如第二轮组合的是:

  • 观察组合的左右两个序数只有中间一位不同,推广到N=2^p都有相同的规律,第 x 轮组合的是右数第 x 位不同的两单元,第 x 位称为第 x 轮的【异位】。

四、组合右边的单元乘以W(k,N),留意到 k 的规律,和组合的两个数密切相关,是两数【异位】右边部分的左移b位,第1轮b=2,第2轮b=1,第3轮b=0。举例:

  • 第二轮,对于组合1-3(二进制:001-011)
  • 异位是右数第二位,异位右边是二进制1
  • 对1左移b=1位,得10,十进制2

五、组合的两个单元执行交叉加减法,左'=左+右,右'=左-右

算法都在上面了,你可以动手写程序了。你需要看清楚图中每一轮转换的规律,写出适应任意N=2^p的计算代码。

【优化】基于大数运算来计算圆周率的算法大部分时间都用来FFT,所以提高FFT性能是关键。下面是一些优化思路,需要一定的编程功底,可根据你自己的能力选择使用哪些。

【优化1:旋转因子缓存】计算中重复使用旋转因子,计算机计算cos和sin的速度较慢,可将其一次计算并保存在数组中,计算时直接提取,在整个计算过程中会被重复调用,可有效节约计算时间。

【优化2:CPU指令集优化】FFT中涉及大量复数运算,SSE系列指令集非常适合对其进行优化。VC和GCC编译器都有提供库(但不统一)来直接调用SSE指令集,经在VC中测试,对32位编译的程序性能提升达3倍,而64位编译程序提升没那么多但仍然显著,可能64位编译时已经有一些智能的方法自动调用SSE。

【优化3:多线程优化】FFT算法是典型的分治算法,非常适合多线程并行运算。注意上面的计算流程图,除了第三轮计算,前面的计算是左右独立的,意味着前两轮可以建立2个线程分别计算左右部分,两者汇合(Join)计算最后一轮。同理,第一轮可以分成4个线程计算,计算完毕后两两汇合成两组计算。当N非常大,可以很容易将计算量分摊到所有CPU线程中,例如N=,为10轮计算,CPU有4个线程,则前8轮分成4个线程计算,第9轮汇合成2个线程计算,最后第10轮一个线程完成计算。

费了九牛二虎之力搞掂了FFT,我们的高速算pi已经几近成功了,接下来我们用牛顿迭代法把除法和开方转成乘法,那就功德完满了。

我们先来复习一下牛顿迭代法:牛顿迭代法可以用来计算任意2阶可导的函数f(x)=0的数值解。
假设方程在附近有一个解,那么使用迭代公式:
那么x序列将越来越接近方程的解,每迭代一次有效小数增加一倍。

对于除法和开方,f(x)有不只一种,要找到适当的f(x),使迭代过程中没有除法和开方,只有乘法。

【除法算法】计算a/b,先用牛顿迭代计算1/b,f(x)是:
迭代式:
【取倒数算法】
迭代结果再乘以a,即为a/b。
由于牛顿跌代法对任意一步的误差都不敏感,你可以取任意初始值来计算,不过最好的方法是将b的前几个单元换算成浮点数,计算1/b然后开始迭代。由于开始迭代时精度低,无需用太高精度的数组,b也截断到较低精度,然后每迭代一次的精度增加一倍,直到达到指定的精度。如此,每迭代一次的开销增加一倍,总计算量的一半在进行最后一轮的迭代。最后一轮迭代有2次乘法,再加上最后计算一个乘法,总计算量约是同样精度乘法的5倍,时间复杂度依然是O(N

【开方算法1】计算a的开方:
迭代(这个方法不好,看下面的【开方算法2】):
 
虽然我们已经能算除法了,但这就成了迭代中的迭代,让人不舒服。

【开方算法2】
先计算:
迭代 x 的f(x)是:
迭代式:
迭代完成后再用上面【取倒数算法】计算1/x即为

计算方法与除法很相似,你甚至可以重用大部分代码。一次迭代3次乘法,开始时只需要截断到低精度计算并逐次精度加倍,总迭代计算量是6个乘法,加上取倒数是4个乘法,总共是10个乘法的计算量,时间复杂度依然是O(N logN)。测试证明【开方2】比【开方1】快20%左右,而且避免了迭代套迭代,推荐使用。

本文主要介绍了这几个计算圆周率相关的要点:

  1. 牛顿迭代高精度除法和开方


事实上如果你掌握了高性能高精度计算方法,可以挑战一些更有难度但更快的算法。

目前世界纪录是用巨型机跑出来的10万亿位,如果用SuperPI和普通CPU计算的话,它要吃掉600T的内存和耗时20年。

笔者使用OpenCL编程开发的程序可使用GPU进行计算,可在3秒内完成一百万位计算,参见:

这个具体可以查一下芯片内核指令代码的手册,如果有硬件乘法指令的系统,16位运算与32位运算所用时间,一般是一样的。

楼上的忽略了32位机这个前提。一般在32位机上是会有硬件乘法指令或者是乘法器部件,是可以在单周期内计算32位乘法的。而16位的乘法,一般都是扩展成32位的乘法来实现,所以说时间一般来说是一样的,都是一个周期。

  1. 没有硬件乘法,使用移位加程序模拟乘法功能,32位的乘法时间可能会比16位乘多一倍,当然这还取决于乘法模拟算法的具体实现。

  2. 用32位实现16位乘法,可能需要对16位数进行额外的扩展调整,这时是16位的慢一点。一般在C语言中就是根据系统字长规定int的类型,所以在C语言中如果不考虑移植,多使用int可能会提高编译的效率,在32位机上最好使用32位的数据进行计算反而更快一些。

本文描述基本的32位X86汇编语言的一个子集,其中涉及汇编语言的最核心部分,包括寄存器结构,数据表示,基本的操作指令(包括数据传送指令、逻辑计算指令、算数运算指令),以及函数的调用规则。个人认为:在理解了本文后,基本可以无障碍地阅读绝大部分标准X86汇编程序。当然,更复杂的指令请参阅Intel相关文档。

主要寄存器如下图所示:

X86处理器中有8个32位的通用寄存器。由于历史的原因,EAX通常用于计算,ECX通常用于循环变量计数。ESP和EBP有专门用途,ESP指示栈指针(用于指示栈顶位置),而EBP则是基址指针(用于指示子程序或函数调用的基址指针)。如图中所示,EAX、EBX、ECX和EDX的前两个高位字节和后两个低位字节可以独立使用,其中两位低字节又被独立分为H和L部分,这样做的原因主要是考虑兼容16位的程序,具体兼容匹配细节请查阅相关文献。

应用寄存器时,其名称大小写是不敏感的,如EAX和eax没有区别。

可以在X86汇编语言中用汇编指令.DATA声明静态数据区(类似于全局变量),数据以单字节、双字节、或双字(4字节)的方式存放,分别用DB,DW, DD指令表示声明内存的长度。在汇编语言中,相邻定义的标签在内存中是连续存放的。

还可以声明连续的数据和数组,声明数组时使用DUP关键字

现代X86处理器具有232字节的寻址空间。在上面的例子中,我们用标签(label)表示内存区域,这些标签在实际汇编时,均被32位的实际地址代替。除了支持这种直接的内存区域描述,X86还提供了一种灵活的内存寻址方式,即利用最多两个32位的寄存器和一个32位的有符号常数相加计算一个内存地址,其中一个寄存器可以左移1、2或3位以表述更大的空间。下面例子是汇编程序中常见的方式

下面是违反规则的例子:

在声明内存大小时,在汇编语言中,一般用DB,DW,DD均可声明的内存空间大小,这种现实声明能够很好地指导汇编器分配内存空间,但是,对于

如果没有特殊的标识,则不确定常数2是单字节、双字节,还是双字。对于这种情况,X86提供了三个指示规则标记,分别为BYTE PTR, WORD PTR, and DWORD PTR,如上面例子写成:mov BYTE PTR [ebx],

汇编指令通常可以分为数据传送指令、逻辑计算指令和控制流指令。本节将讲述其中最重要的指令,以下标记分别表示寄存器、内存和常数。

任何8位、16位或32位常数

mov指令将第二个操作数(可以是寄存器的内容、内存中的内容或值)复制到第一个操作数(寄存器或内存)。mov不能用于直接从内存复制到内存,其语法如下所示:

push指令将操作数压入内存的栈中,栈是程序设计中一种非常重要的数据结构,其主要用于函数调用过程中,其中ESP只是栈顶。在压栈前,首先将ESP值减4(X86栈增长方向与内存地址编号增长方向相反),然后将操作数内容压入ESP指示的位置。其语法如下所示:

pop指令与push指令相反,它执行的是出栈的工作。它首先将ESP指示的地址中的内容出栈,然后将ESP值加4. 其语法如下所示:

 lea实际上是一个载入有效地址指令,将第二个操作数表示的地址载入到第一个操作数(寄存器)中。其语法如下所示:

3.2 算术和逻辑指令

add指令将两个操作数相加,且将相加后的结果保存到第一个操作数中。其语法如下所示:

sub指令指示第一个操作数减去第二个操作数,并将相减后的值保存在第一个操作数,其语法如下所示:

inc,dec分别表示将操作数自加1,自减1,其语法如下所示:

整数相乘指令,它有两种指令格式,一种为两个操作数,将两个操作数的值相乘,并将结果保存在第一个操作数中,第一个操作数必须为寄存器;第二种格式为三个操作数,其语义为:将第二个和第三个操作数相乘,并将结果保存在第一个操作数中,第一个操作数必须为寄存器。其语法如下所示:

idiv指令完成整数除法操作,idiv只有一个操作数,此操作数为除数,而被除数则为EDX:EAX中的内容(一个64位的整数),操作的结果有两部分:商和余数,其中商放在eax寄存器中,而余数则放在edx寄存器中。其语法如下所示:

逻辑与、逻辑或、逻辑异或操作指令,用于操作数的位操作,操作结果放在第一个操作数中。其语法如下所示:

位翻转指令,将操作数中的每一位翻转,即0->1, 1->0。其语法如下所示:

位移指令,有两个操作数,第一个操作数表示被操作数,第二个操作数指示位移的数量。其语法如下所示:

X86处理器维持着一个指示当前执行指令的指令指针(IP),当一条指令执行后,此指针自动指向下一条指令。IP寄存器不能直接操作,但是可以用控制流指令更新。

一般用标签(label)指示程序中的指令地址,在X86汇编代码中,可以在任何指令前加入标签。如:

如第二条指令用begin指示,这种标签的方法在某种程度上简化了汇编程序设计,控制流指令通过标签实现程序指令跳转。

控制转移到label所指示的地址,(从label中取出执行执行),如下所示:

条件转移指令,条件转移指令依据机器状态字中的一些列条件状态转移。机器状态字中包括指示最后一个算数运算结果是否为0,运算结果是否为负数等。机器状态字具体解释请见微机原理、计算机组成等课程。语法如下所示:

cmp指令比较两个操作数的值,并根据比较结果设置机器状态字中的条件码。此指令与sub指令类似,但是cmp不用将计算结果保存在操作数中。其语法如下所示:

比较var指示的4字节内容是否为10,如果不是,则继续执行下一条指令,否则,跳转到loop指示的指令开始执行

这两条指令实现子程序(过程、函数等意思)的调用及返回。call指令首先将当前执行指令地址入栈,然后无条件转移到由标签指示的指令。与其它简单的跳转指令不同,call指令保存调用之前的地址信息(当call指令结束后,返回到调用之前的地址)。

ret指令实现子程序的返回机制,ret指令弹出栈中保存的指令地址,然后无条件转移到保存的指令地址执行。

call,ret是函数调用中最关键的两条指令。具体细节见下面一部分的讲解。语法为:

为了加强程序员之间的协作及简化程序开发进程,设定一个函数调用规则非常必要,函数调用规则规定函数调用及返回的规则,只要遵照这种规则写的程序均可以正确执行,从而程序员不必关心诸如参数如何传递等问题;另一方面,在汇编语言中可以调用符合这种规则的高级语言所写的函数,从而将汇编语言程序与高级语言程序有机结合在一起。

调用规则分为两个方面,及调用者规则和被调用者规则,如一个函数A调用一个函数B,则A被称为调用者(Caller),B被称为被调用者(Callee)。

下图显示一个调用过程中的内存中的栈布局:

在X86中,栈增长方向与内存编号增长方向相反。

调用者规则包括一系列操作,描述如下:

1)在调用子程序之前,调用者应该保存一系列被设计为调用者保存的寄存器的值。调用者保存寄存器有eax,ecx,edx。由于被调用的子程序会修改这些寄存器,所以为了在调用子程序完成之后能正确执行,调用者必须在调用子程序之前将这些寄存器的值入栈。

2)在调用子程序之前,将参数入栈。参数入栈的顺序应该是从最后一个参数开始,如上图中parameter3先入栈。

3)利用call指令调用子程序。这条指令将返回地址放置在参数的上面,并进入子程序的指令执行。(子程序的执行将按照被调用者的规则执行)

当子程序返回时,调用者期望找到子程序保存在eax中的返回地址。为了恢复调用子程序执行之前的状态,调用者应该执行以下操作:

2)将栈中保存的eax值、ecx值以及edx值出栈,恢复eax、ecx、edx的值(当然,如果其它寄存器在调用之前需要保存,也需要完成类似入栈和出栈操作)

如下代码展示了一个调用子程序的调用者应该执行的操作。此汇编程序调用一个具有三个参数的函数_myFunc,其中第一个参数为eax,第二个参数为常数216,第三个参数为var指示的内存中的值。

在调用返回时,调用者必须清除栈中的相应内容,在上例中,参数占有12个字节,为了消除这些参数,只需将ESP加12即可。

 _myFunc的值保存在eax中,ecx和edx中的值也许已经被改变,调用者还必须在调用之前保存在栈中,并在调用结束之后,出栈恢复ecx和edx的值。

被调用者应该遵循如下规则:

1)将ebp入栈,并将esp中的值拷贝到ebp中,其汇编代码如下:

上述代码的目的是保存调用子程序之前的基址指针,基址指针用于寻找栈上的参数和局部变量。当一个子程序开始执行时,基址指针保存栈指针指示子程序的执行。为了在子程序完成之后调用者能正确定位调用者的参数和局部变量,ebp的值需要返回。

2)在栈上为局部变量分配空间。

4)在上述三个步骤完成之后,子程序开始执行,当子程序返回时,必须完成如下工作:

  4.1)将返回的执行结果保存在eax中

  4.4)通过弹出栈中保存的ebp值恢复调用者的基址寄存器值。

  4.5)执行ret指令返回到调用者程序。

子程序首先通过入栈的手段保存ebp,分配局部变量,保存寄存器的值。

在子程序体中,参数和局部变量均是通过ebp进行计算。由于参数传递在子程序被调用之前,所以参数总是在ebp指示的地址的下方(在栈中),因此,上例中的第一个参数的地址是ebp+8,第二个参数的地址是ebp+12,第三个参数的地址是ebp+16;而局部变量在ebp指示的地址的上方,所有第一个局部变量的地址是ebp-4,而第二个这是ebp-8.

参考资料

 

随机推荐