算法系列：求幂算法

算法系列：求幂算法

2014/08/06 | 分类： IT技术 | 1 条评论 | 标签： TAOCP, 算法, 编程

本文由 伯乐在线 - 知曰 翻译自 Alvaro Videla。未经许可，禁止转载！
欢迎加入：技术翻译小组，或分享原创到伯乐头条。

1.快速求幂算法

在这篇文章我会展示怎样通过求一个数的幂的基本思路，来引导我们发现一些抽象的东西比如半群和含幺半群。

有一个很有名的对一个数求幂的算法,也就是说,求一个数x的n次方或者这样简单表示：x^n。Donald Knuth在TAOCP的4.63节 求幂值中提出这个算法。

这个算法很简单的实现就是x乘以自己n次，但是在这里当然会提供一种比这种方式更快的算法。正在谈论的算法通常被称作二进制法(binary method)、梯度求幂(the powering ladder)或者反复平方法(repeated-squaring algorithm)

假设我们想计算2^23，在这里x = 2，n = 23，这个算法首先把23表示成二进制的形式10111。扫描这个二进制数(10111)每当遇到0或1，则相应的求x的平方或者乘以x。

这个方法有一个问题就是它扫描二进制表示的数是从左到右进行的，但是对于计算机通常以相反的方向能够更容易实现，因此Knuth提出一个替代的算法。

一个出自TAOCP的4.63节的算法A的简单实现如下：

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

function power1( $x , $n ) {      $y = 1;        while (true) {          $t = $n % 2;          $n = floor ( $n /2);            if ( $t == 1) {              $y = $y * $x ;          }            if ( $n == 0) {              break ;          }            $x = $x * $x ;      }        return $y ; }

这个函数需要两个整数，$x和$n然后返回$x的$n次幂作为结果。

首先创建一个辅助变量$y并且初始化为1，把它作为乘法的主体。

然后函数在每次循环迭代的时候扫描$n的二进制表示的数。如果遇到1则$y乘上$x，然后赋值回$y。每次循环都会计算$x的平方，并且把它赋值回$x。

遇到1意味着当前$n的值不能被2整除，换句话说就是，$n % 2 == 1。

同样的每次循环$n都会折半，然后向下取整得到结果。当$n等于0的时候，我们结束循环并且返回$y的值。

函数power能够这样被调用:

1 2	1024 == power1(2, 10); => true

我能想象你现在就像这个gif中的男孩。

尽管这个算法看起来像一个 “呵呵，真有意思” (无语了，你别说了，我根本不关心)的故事，实际上当它用来计算非常的大数时时十分高效的。例如有很多的素数测试算法都是依赖这个算法的不同变式。

2.增加一些抽象

到目前为止还没有什么意想不到的事情发生，但是如果我们注意到求一个数的幂实际上和一个数自乘多次是等价的，我们也可以看到乘法实际上等价于自加多次。举个例子2 * 5能够像这样被计算2 + 2 + 2 + 2 + 2。

我们能把这个算法转换成一种更普遍的形式使它能同样应用在乘法还有加法上吗?当然可以，我们仅仅需要改变几样东西。

在当前实现中，我们创建$y作为乘法的主体，并设置为1。如果我们想把算法用在加法上，我们需要把$y设置为0。因此我们仅需要改变函数的单位元素的值。

第二步要提供一个函数给我们的算法，它能够作乘法或者加法。为了实现这个目的我们会传递一个担当二元运算的函数。例如：一个需要两个参数的函数。这个函数需要遵循以下的规则。必须满足：a·（ b · c ） = (a · b ) · c。还要求返回结果的类型必须和两个输入参数的类型一致。

幸运的是加法和乘法都满足结合律，因此我们能够仅在一个函数中包含他们然后把它传递给我们的power算法。

这里是这个算法新的实现：

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

function power2( $x , $n , $id , $f ) {        $y = $id ;          while (true) {            $t = $n % 2;            $n = floor ( $n /2);              if ( $t == 1) {                $y = $f ( $y , $x );            }              if ( $n == 0) {                break ;            }              $x = $f ( $x , $x );        }          return $y ;    }

我们能够像这样调用它:

1 2	1024 == power2(2, 10, 1, function ( $a , $b ) { return $a * $b ; });   => true

记住传递进我们算法的运算必须是可结合的,举个例子,减法不能被用在这里由于10 - ( 5 - 3) = 8但是(10 - 5 ) - 3 = 2。

3.附加更抽象的概念

从数学的角度说这个算法能够在任何满足结合律的代数结构中有效(在这个案例中就是整数的乘法和加法)，换言之，它能够用在半群中，引用一本关于群论的书。

一个半群的集合S含有一个可结合的运算 · ；
也就是说,x·(y · z) = (x · y) · z 对于所有的x, y, z ∈ S都成立。

同样，这个集合必须有一个单位元素使得它有一个独异点:

一个独异点是一个集合M含有一个可结合运算·；伴有一个单位元素e∈ M满足e·x = x· e =  x对于所有x∈ M都成立。

在这个预设条件下，有什么我们经常用在编程上的结构能使用这个算法的呢？如果你是一个web开发者，你不需要费大力气去获取strings。对于字符串（strings），使用string append作为二元操作而且空字符串(empty string)作为单位元素同样会带来类似的结果。如果一个字符串想重复n次,我们创建下面的函数:

1 2 3 4 5

function repeat( $s , $n ) {      return power2( $s , $n , "" , function ( $a , $b ) {                 return $a . $b ;             }); }

测试:

1 2	"aaaaaaaaaa" == repeat( "a" , 10); => true

现在考虑一下数组(arrays)(或者其它语言称为列表(lists))。我们想把一个数组复制n次。在这里空数组是单位元素，对PHP来说array_merge会用来作为二元操作。

1 2 3 4 5

function repeat_el( $el , $n ) {      return power2( array ( $el ), $n , array (), function ( $a , $b ) {              return array_merge ( $a , $b );             }); }

结果:

1 2 3

$arr = repeat_el( "a" , 10); 10 == count ( $arr ); => true

从上不难看出,像求一个数幂运算的这样简单事情给我们带来一个优雅的算法,它能被运用一些事情上,像重复的东西还有数组里的元素。

4.延伸阅读

• 这里的快速求幂算法是基于TAOCP中，卷二的4.63节。
• 所有的关于工作原理的解答都可以在TAOCP或者在这本书《A Computational Introduction to Number Theory and Algebra》上找到,这本书的PDF版本在作者的主页上可以免费下载。浏览章节：“Computing with large integers - The repeated squaring algorithm”
• 如果你想学习这个算法的一些用法或者想知道更多这个算法背后的理论,请查阅这本叫做《Elements of Programming》的书。这本书非常了不起，它定义了不同类型的函数和使用类型系统确定函数是否是可结合的，二元的等等。作者是C++STL的设计者,所以这本书的内容可能会比较理论化,然后它能够直接应用在面向对象编程(OOP)。
• 半群 和 含幺半群的引用来自于《Handbook of Computational Group Theory.》。一本非常有趣的书，如果你对计算群论有兴趣的话。
• 如果你想学习更多有关幺半群还有它们的实现。《Learn You a Haskell》里的有个章节非常有趣的介绍它:Functors, Applicative Functors and Monoids
• 这是一个十分有趣的练习,通过实现这些概念使用PHP和OOP,对于不喜欢使用PHP无爱的人,也可以选择其它你喜欢的语言。

5.你是想说Haskell?

既然我已经提及一本Haskell的书，这里有一个Haskell实现的求幂算法，使用的递归算法来自于这本书《Prime Numbers: A Computational Perspective》

1 2 3 4 5 6 7

power :: (Eq a, Integral b) => (a -> a -> a) -> a -> b -> a power f a n    | n == 1 = a    | even n = square a (n `div` 2)    | otherwise = f a (square a ((n-1) `div` 2))    where      square a ' n' = f (power f a ' n' ) (power f a ' n' )

几个函数调用的结果：

1 2 3 4 5 6 7 8 9

*Main> :load pow.hs    [1 of 1] Compiling Main             ( pow.hs, interpreted )    Ok, modules loaded: Main.    *Main> power (*) 2 10    1024    *Main> power (+) 2 10    20    *Main> power (++) "a" 10    "aaaaaaaaaa"

正如你所看到的,这个函数调用一个function(a->a->a)，例子中，对于integers使用*或者+,对于lists使用++。

我希望你会觉得这边文章有趣或者激起你学习与编程有关的数学的欲望。因为我认为我们掌握得越多数学方面的知识，我们就能更好的使用抽象的东西。