Chapter 02 The First Algorithm

摩西快速学习了算术和几何……
这些知识他是从埃及人那里学到的。埃及人研习数学和其它学问。

《形而上学》亚里士多德

算法 是用以完成某项计算任务的有限步骤。算法跟计算机编程有紧密的联系,以至于大多数知道这个词的人以为算法的使用是从计算机科学开始。但实际上,人们使用算法已经有几千年的历史。数学中充满了算法,有一些我们到今天还在使用。学生们学习的长加法就是一种算法。

尽管使用算法的历史很长,算法这个概念并不是一开始就存在,它是后来发明出来的。尽管我们不知道是谁首先发明了这个概念,却明确知道有些算法至迟在4000多年前的古埃及就存在了。

 

古埃及文明的核心是尼罗河,因为其农业有赖于尼罗河泛滥带来的肥沃泥土。问题是,尼罗河一发水,用来划分权属的界标也都被冲掉了。埃及人惯用绳子测量长度,因此发展出一套方法,可以根据书面记录重建界标。一群制定的祭司专司此事,他们学过这些数学技术,因此成为“拉绳子的人(rope-strecher)”。希腊人后来会把他们叫做几何学家(geometers),意思是“土地丈量者”。

很不幸我们并没有多少古埃及数学知识的文字记录,当时的数学文件仅有两篇留存至今。我们关心的一个,叫做Rhind Mathematical Papyrus,名字来自于19世纪在埃及买到它的苏格兰收藏者。这篇文献成于约公元前1650年,抄写者名叫Ahmes。文献中有一系列算术和几何问题,还有一些辅助计算的表格。其中包含一个快速乘法技术和一个快速除法技术,是最早的有记录的算法。我们首先来看一看这个快速乘法算法,(我们很快会看到)它至今仍然是重要的计算技术。

2.1 埃及乘法

跟其它古文明一样,埃及人使用的数字系统中没有数位概念,也没有零。因此,乘法非常复杂,只有少数专家能够掌握。(想象一下只使用罗马数字进行两位数乘法)

首先,什么是乘法?粗略地说,乘法就是“把一个数跟它自身反复相加”。严格一些,我们可以首先把乘法分为两种情况:乘以1,或者乘以一个大于1的数。

乘数为1的乘法的定义是:

$$1a=a$$                                   (2.1)

然后,我们考虑如何计算比我们已经计算出来的再多一倍的情形。有些读者可以意识到这是一个递推,我们会在后面正式使用这个技术。

$$(n+1)a=na+a$$                 (2.2)

计算na的厂家做法是把n个a加起来。但是,在n很大的时候,这可能非常繁琐,因为需要进行n-1次加法。使用C++,这个算法是这样的:

int multiply0(int n, int a){
    if (n==1) return a;
    return multiply0(n-1, a)+a;
}

这样两行代码对应于(2.1)和(2.2)。a和n都必须是正数,古埃及当时的约定也是如此。

Ahmes所描述的算法,这个算法在古希腊被称为“埃及乘法”,而很多现代作者则称之为“俄罗斯农民算法”,建立在如下洞察上:

$$ 4a=((a+a)+a)+a=(a+a)+(a+a) $$

这个优化利用了加法的结合律:

$$a+(b+c)=(a+b)+c $$

这样,我们只需要计算一次a+a,减少了计算次数。

核心思想在于,不断地把n取半,同时把a加倍,使用数的2的次方倍相加来计算乘法结果。在当时,人们不会使用变量a和n来描述算法,而是由作者给出示例,然后说“其它数字同理”。Ahmes也是这么做的,他使用下面的表格来演示这种算法,n=41,a=59:

1        ♦         59
2                  118
4                  236
8        ♦        472
16                944
32      ♦      1888

每个条目的最左边是2的幂,右边是前一条目的数字的倍数(数字加倍相对容易计算)。有标记的值对应41的二进制表示中数值为1的位。这个表格的意思是说:

$$41*59=(1*59)+(8*59)+(32*59)$$

其中右边的每一个乘积可以由对59加倍响应的次数而得到。

这个算法需要检查n是奇数还是偶数,所以我们可以推断古埃及人知道它们的区别,尽管我们没有直接的证据,但是古希腊人,他们自称从古埃及人那里学到数学,显然知道奇偶的区别。下面是他们的定义方法用现代数学标记表达出来:

$$n=n/2+n/2$$    说明n是偶数
$$n=(n-1)/2+(n-1)/2+1$$  说明n是奇数

我们也将使用这个特点:

odd(n) 意味着 half(n)=half(n-1)

下面是埃及乘法的C++实现:

int multiply1(int n, int a){
    if (n==1) return a;
    int result=multiply1(half(n), a+a);
    if (odd(n)) result=result +a;
    return result;
}

odd(x)很容易实现,测试x的最低位就可以了,half(x)则可以通过对x进行一次右移实现:

bool odd(int n) { return n&0x1;}
int half(int n) { return n>>1; }

multiply1需要进行多少次加法运算?每次调用这个函数,我们需要在a+a的地方进行一次加法,因为我们是不断对n取半,我们需要调用这个函数logn次。在某些时候,我们还需要在result+a的地方进行加法运算,因此加法运算的总次数为:

$$\#(n)=logn + v(n)-1$$

其中v(n)是n的二进制表达中1的个数(population count或者pop count)。由此我们把一个O(n)算法优化成了O(logn)算法。

这个算法已经最优了吗?不见得。比如,如果我们要乘以15,使用之前的公式,我们需要

$$\#(15)=3+4-1=6$$

进行6次加法。但是我们使用下面的办法,就只需要5次加法运算:

int multiply_by_15(int a){
    int b=(a+a)+a;    //b == 3*a
    int c+b+b;        //c == 6*a
    return (c+c)+b;   //12*a + 3*a
}

这样一系列加法运算又称为加法链。这里我们发现了15的最优加法链,尽管如此,Ahmes的算法在多数情况下已经足够好了。

练习2.1

为n<100的数找到最优加法链

有读者可能已经注意到其中有些乘法运算可以更快,只要我们在第一个参数大于第二个参数的情况下交换一下参数顺序就可以了(比如,3*15会比15*3更容易计算)。的确如此,而且古埃及人也知道这一点。但是我们不会在这里进行这样的优化,因为到第七章我们会发现,我们最终要把算法推广到不同的数据类型。在这种情况下,参数的顺序是不能交换的。

2.2 改进算法

函数multiply1的加法运算次数还可以,但是它还要进行logn次递归调用。鉴于函数调用的开销很高,我们希望能够优化程序以避免这个开销。

我们将要利用的一个原则是:It is often easier to do more work rather than less. 特别地,我们准备计算

r+na

其中r是加法计算的中间结果,也就是说,我们要实现一个乘法累计函数,而不仅仅是乘法函数。这种思维方式不仅仅适用于程序编写,也适用于硬件设计和数学——比如在数学中,时常有一般性判断比特定判断更容易证明的情形。

这里是我们的乘法累计函数:

int mult_acc0{int r, int n, int a) {
    if (n==1) return r+a;
    if (odd(n)){
        return mult_acc0(r+a, half(n), a+a);
    }else{
        return mult_acc0(r, half(n), a+a);
    }
}

这里有一个不变形式:r+na=r0+n0a0,其中r0、n0和a0是r、n和a这些变量的初始值。

我们可以通过简化递归获得进一步优化。注意在上面的两行递归调用代码只有第一个参数不同,这样对于奇数情况和偶数情况分别进行递归,不如我们在进行递归前根据情况修改第一个参数的值:

int mult_acc1(int r, int n, int a) {
    if (n==1) return r+a;
    if (odd(n)) r=r+a;
    return mult_acc1(r, half(n), a+a);
}

现在我们的函数是尾递归函数——递归仅仅发生在返回值里面——我们马上会利用这一特点。

我们注意到两点:

  • n=1的情况很少发生;
  • 如果n是偶数,则完全没有必要判断它还是不是1.

因此,我们可以首先判断n的奇偶性,这样把n与1进行比较的次数降低一半:

int mult_acc2(int r, int n, int a) {
    if (odd(n)) {
        r=r+a;
        if (n==1) return r;
    }
    return mult_acc2(r, half(n), a+a);
}

有些程序员认为编译器可以为我们做这样的优化,其实不太可能。编译器不能把一种算法优化成另外一种算法。

到目前为止都很好,不过我们的目标是避免递归以节省函数调用的开销。要达到目的,如果我们的函数是完全尾递归就好了。

定义2.1 一个完全尾递归函数是一个所有递归调用的形式参数都跟函数本身一致的函数。

同样地,我们可以通过在递归调用之前给变量重新赋值来达到这个要求:

int mult_acc3 (int r, int n, int a) {
    if (odd(n)) {
        r=r+a;
        if (n==1) return r;
    }
    n=half(n);
    a=a+a;
    return mult_acc3(r,n,a);
}

到了这里,把它转化成一个循环就很简单了,只需要使用一个while(true)结构替换掉尾递归就可以了:

int mult_acc4(int r, int n, int a) {
    while (true) {
        if (odd(n)) {
            r=r+a;
            if (n==1) return r;
        }
        n=half(n);
        a=a+a;
    }
}

使用这个全新的乘法累加函数,我们可以写一个新版的乘法函数:

int multiply2(int n, int a) {
    if (n==1) return a;
    return mult_acc4(a, n-1, a);
}

注意我们在调用mul_acc4的时候预先把累加结果置为a,这样就减少了一次循环。

效果很不错,除了在n是2的幂的时候。由于我们第一次调用就把n减去1,这样mult_acc4函数被调用的时候第二个参数是一个二进制表达是全1的数,这对于我们的算法是最坏情况。因此,我们要在n是偶数的情况下做一些预处理,以避免这种情形。我们把n取半(同时把n加倍),直到n变成奇数。

int multiply3(int n, int a) {
    while (!odd(n)) {
        a=a+a;
        n=half(n);
    }
    if (n==1) return a;
    return mult_acc4(a, n-1, a);
}

但这样的话,mul_acc4就没有必要再检查n=1的情况了,因为n是奇数,n-1是偶数。因此,我们可以在调用的时候直接把n-1取半,把a再加倍,这样就是我们的最终版:

int multiply4(int n, int a) {
    while (!odd(n)) {
        a=a+a;
        n=half(n);
    }
    if (n==1) return n;
    return mult_acc4(a, half(n-1), a+a);
}

重写代码

在不断变换我们的乘法算法的过程中,可以看到重写代码是很重要的。没有谁一上来就能写出很好的代码,很多时候需要反复尝试才能找到最高效或者最通用的方法。程序员不应该有一遍就过的想法。

在重写过程中你可能会想,一条指令的多少不见得能造成多大的区别。但是也可能你的代码在后来很长时间内被重用很多次(实际上,临时处理措施经常变成代码里面存活时间最长的部分),而且,本来可以省掉的操作现在不起眼,留下来也可能会被修改成非常开销非常高的代码。

努力追求效率的另外一个好处是,在此过程中你不得不逼迫自己深入理解问题本身。与此同时,理解深入必然带来更有效率的实现,这是一个良性循环。

2.3 本章的思考

初等代数的学生都会学习如何通过变换表达式的形式来简化它们。在我们实现埃及乘法的不同实现中,我们也经历了类似的变换过程,重新安排代码语句,使之更清晰更高效。每个程序员都应该养成尝试对代码进行变换的习惯,直到得到最佳形式。

本章我们看到了数学在古埃及是怎样出现的,如何留下了第一个有记录的算法。我们要到靠后面的内容才会回到这个算法并对其进行推广。接下来,我们要向后跳跃1000年,领略一些古希腊人的数学发现。