银河

SKYIV STUDIO

  博客园 :: 首页 :: 博问 :: 闪存 :: :: :: 订阅 订阅 :: 管理 ::
  268 随笔 :: 2 文章 :: 2606 评论 :: 48 引用

公告

最近一段时间,我在 Timus Online Judge 网站做 ACM 题。发现其中不少题目需要用到 BigInteger,例如:

于是,我就开始寻找可以用于 C# 的 BigInteger 类。结果如下表所示:

项目 速度 自动增长 进制 ModPowPowSqrt 源文件大小 备注
Fcl35Biginteger 最快 232 1,755行,47.1 KB 是 struct,不是 class
SkyivBigInteger
109 335行,9.16 KB 源文件短小精悍,功能全
ChewBigInteger
232 3,335行,102 KB 源文件中有大量的注释
JavaBigInteger 极慢 28? (vjslib.dll)
需要 J# run-time library

在上表中,JavaBigInteger 是指 java.math.BigInteger,她在 vjslib.dll 中,需要 J# run-time library (Visual Studio 2005/2008 中已经包括 J# run-time library)。用 .NET Reflector 打开 C:\Windows\Microsoft.NET\Framework\v2.0.50727\vjslib.dll 文件,找到 java.math.BigInteger,查看其源程序代码如下:

1 [Serializable, JavaInterfaces("2;System/IConvertible;java/lang/Comparable;")]
2 public class BigInteger : Number, IConvertible, Comparable
3 {
4   // Fields
5   private BigInteger __absOfthis;
6   [JavaFlags(0x1002)]
7   internal sbyte[] __bits;
8   // 
9 }

根据第 7 行的 internal sbyte[] __bits; 猜想 java.math.BigInteger 的内部实现是基于 28 进制的。


ChewBigInteger 就是 Chew Keong TANCodeProject 中的 C# BigInteger Class 。这是个很有名的 C# Biginteger 类,它写于2002年8月到9月,已经有很久的历史了。但是她有个致命的缺点,就是不能自动增长,需要事先指定最大允许有多大的数,如果在运算过程中产生超过你事先指定大小的数,就会引发异常。

她的部分源程序代码如下所示:

 1 namespace ChewKeongTAN
 2 {
 3   /*public */class BigInteger
 4   {
 5     // maximum length of the BigInteger in uint (4 bytes)
 6     // change this to suit the required level of precision.
 7 
 8     private const int maxLength = 26591// == Math.Ceiling ( digits / log10(2^32) )
 9 
10     private uint[] data = null;          // stores bytes from the Big Integer
11     public int dataLength;               // number of actual chars used
12     
13     // 
14     
15   }
16 }

她内部使用一个 uint[] 来存储,是基于 232 ≈ 4.3 x 109 进制的。需要事先设定 maxLength 的值(程序中第 8 行),用以指定最多允许多少个 232 进制的数字,换算用的公式是:

maxLength = Math.Ceiling ( digits / log10(232) ) ≈ digits / 9.6

上式中 digits 表示十进制数字的个数。


Fcl35Biginteger 就是 .NET Framework 3.5 Base Class Library 中的 System.Numeric.BigInteger,她在 System.Core.dll 中。注意,她是一个 internal struct。首先,她是一个 struct,不是一个 class,也就是说她是值类型而不是引用类型。其次,致命的是,她是 internal 的,而不是 public 的,也就是说,她就在那里,但是你就是不能使用她。实际上,在 .NET Framework 3.5 的 Beta 和 latest CTP 版本中,她还是 public 的,只不过是在 .NET Framework 3.5 正式发布时变成了 internal 的。这里有一篇 BCL Team 的 Melitta Andersen 于2008年1月4日发表的 Blog: Where did BigInteger go? 谈到这件事,要点是:

So why was BigInteger cut?  The basic rationale behind making BigInteger internal was that it just wasn't ready to ship.  We thought our implementation met the needs for a BigInteger type.  But then we had some other teams take a look at it and they pointed out some performance and compatibility issues that we just didn't have time to fix before we shipped.

于是, .NET Reflector 又上场了,用她打开 C:\Windows\assembly\GAC_MSIL\System.Core\3.5.0.0__b77a5c561934e089\System.Core.dll 文件,找到 System.Numeric.BigInteger,在 Disassembler 窗口中,点击最下方的 Expand Methods,然后把源程序代码全部复制到 Visual Studio 2008 的 IDE 窗口(注意,目标 Framework 设定为 .NET Framework 2.0)中。按下编译按钮,天啊,出现了好多错误。于是,就开始修改源程序代码:

  • 实现两个空的 sealed class: ImmutableAttribute 和 PureAttribute。(前一个好理解,后一个不知是什么东东。当然,实际上她们不应该是空的)
  • 实现一个 static class: Contract,其中包括两个重载的 public static void Requires 方法。(当然,实际上的 Contract.Requires 应该更复杂)
  • 实现一个 static class: Res,其中包括若干个 public static readonly string 字段。(当然,实际上 Res 应该从资源文件中获得本地化的数据)
  • 将程序中出现的某些常数值还原为 private const 字段。注意,这不是必须的,只是为了使程序更容易看懂。如果不做这一步,就可以将程序中所有的 private const 字段全部删除。
  • 修改程序中某些变量的名称。注意,这也不是必须的,也只是为了使程序更容易看懂。

然后,再次编译,耶,通过了。但是,还是高兴得太早了一点,运行一些测试,发现答案根本不对。

经检查,发现是 .NET Reflector 在反汇编为 C# 代码时会丢失一些强制类型转换,这在我的上一篇随笔 浅谈 Math.BigMul 方法 中已经提到过了:她把 Math.BigMul 方法中的唯一的语句 return (long)a * b; 错误地反汇编为 return (a * b); 了,但是反汇编为 IL 代码是正确的。

于是,只好逐条语句检查,如果需要进行强制类型转换,就自己加上,同时还去掉一些不必要的强制类型转换和括号,以便使源程序代码更易读。真是一个很烦人和困难和工作,热切盼望 .NET Reflector 能够尽快修复这个 BUG。:)

好了,经过调试,程序终于运行正常了。下面就是修改后的部分源程序代码:

 1 namespace Fcl35.Numeric
 2 {
 3   [Serializable, StructLayout(LayoutKind.Sequential), Immutable, ComVisible(false)]
 4   public /*internal*/ struct BigInteger : IFormattable,
 5     IEquatable<BigInteger>, IComparable<BigInteger>, IComparable
 6   {
 7     private const ulong Base = uint.MaxValue + 1UL// 0x100000000L;
 8     private readonly short _sign;
 9     private readonly uint[] _data;
10     private int _length;
11     // 
12   }
13 }

ChewKeongTAN.BigInteger 类一样,Fcl35.Numberic.BigInteger 结构内部也使用一个 uint[] 来存储,同样也是基于 232 ≈ 4.3 x 109 进制的。但是,最重要的一点,Fcl35.Numberic.BigInteger 结构是动态增长的,不需要事先指定最大允许有多大的数。


现在来讨论我最近写的一个 Skyiv.Numeric.BigInteger 类,她是从我在做 ACM 题时写的极简单的如下所示的 BigInteger 类发展而来的:

1 sealed class BigInteger
2 {
3   int[] digits = new int[1800];
4   // 
5 }

这个极简单的 BigInteger 类的全部源程序代码可以在本随笔开头给出的 URL 中找到,只有五十多行。她是基于 10 进制的,内部使用一个 int[] 来存储,需要事先指定该数组的大小,不能动态增长,而且只能表示非负整数。 

下面是 Skyiv.Numeric.BigInteger 的源程序代码,非常短小精悍,只有区区 335 行,该有的功能基本上都有了。实现了基本的算术运算和逻辑运算,还实现了 Pow 和 Sqrt 方法。

  1 using System;
  2 using System.Text;
  3 using System.Collections.Generic;
  4 
  5 namespace Skyiv.Numeric
  6 {
  7   sealed class BigInteger : IEquatable<BigInteger>, IComparable<BigInteger>
  8   {
  9     static readonly int Len = 9// int.MaxValue = 2,147,483,647
 10     static readonly int Base = (int)Math.Pow(10, Len);
 11 
 12     int sign;
 13     List<int> data;
 14 

她内部使用一个 List<int> 来存储(程序中第 13 行),是动态增长的,基 于 109 进制(第 9 和 10 行)。

 15     BigInteger(long x)
 16     {
 17       sign = (x == 0? 0 : ((x > 0? 1 : -1);
 18       data = new List<int>();
 19       for (ulong z = Abs(x); z != 0; z /= (ulong)Base) data.Add((int)(z % (ulong)Base));
 20     }
 21 
 22     BigInteger(BigInteger x)
 23     {
 24       sign = x.sign;  // x != null
 25       data = new List<int>(x.data);
 26     }
 27 
 28     public static implicit operator BigInteger(long x)
 29     {
 30       return new BigInteger(x);
 31     }
 32 

她只有两个构造函数,而且都是 private 的。要获得该类的实例,主要途径之一是通过从 long、int 等类型隐式转换(第 28 到 31 行)而来,例如 BigInteger a = 2; 语句。

 33     public static BigInteger Parse(string s)
 34     {
 35       if (s == nullreturn null;
 36       s = s.Trim().Replace(","null);
 37       if (s.Length == 0return 0;
 38       BigInteger z = 0;
 39       z.sign = (s[0== '-'? -1 : 1;
 40       if (s[0== '-' || s[0== '+') s = s.Substring(1);
 41       int r = s.Length % Len;
 42       z.data = new List<int>(new int[s.Length / Len + ((r != 0? 1 : 0)]);
 43       int i = z.data.Count - 1;
 44       if (r != 0) z.data[i--= int.Parse(s.Substring(0, r));
 45       for (; i >= 0; i--, r += Len) z.data[i] = int.Parse(s.Substring(r, Len));
 46       z.Shrink();
 47       return z;
 48     }
 49 

获得该类的实例的主要途径之二是调用静态的 Parse 方法,例如:var a = BigInteger.Parse("-123,456,789"); 语句。

 50     public static void Swap(ref BigInteger x, ref BigInteger y)
 51     {
 52       BigInteger z = x;
 53       x = y;
 54       y = z;
 55     }
 56 
 57     public static ulong Abs(long x)
 58     {
 59       return (x < 0? (ulong)-x : (ulong)x;
 60     }
 61 
 62     public static BigInteger Abs(BigInteger x)
 63     {
 64       if (x == nullreturn null;
 65       BigInteger z = new BigInteger(x);
 66       z.sign = Math.Abs(x.sign);
 67       return z;
 68     }
 69 

上面三个静态的公共方法都很简单,主要是为了方便,如果不提供这些方法也没关系。要注意的是第二个 Abs(long) 方法返回 ulong 类型,而且从不抛出异常。而 Math.Abs(long) 返回 long 类型,而且当参数为 long.MinValue 时会抛出异常。

 70     public static BigInteger Pow(BigInteger x, int y)
 71     {
 72       if (x == nullreturn null;
 73       BigInteger z = 1, n = x;
 74       for (; y > 0; y >>= 1, n *= n) if ((y & 1!= 0) z *= n;
 75       return z;
 76     }
 77 
 78     public static BigInteger operator +(BigInteger x)
 79     {
 80       if (x == nullreturn null;
 81       return new BigInteger(x);
 82     }
 83 
 84     public static BigInteger operator -(BigInteger x)
 85     {
 86       if (x == nullreturn null;
 87       BigInteger z = new BigInteger(x);
 88       z.sign = -x.sign;
 89       return z;
 90     }
 91 

注意,重载的单目 + 操作符不仅是为了和重载的单目 - 操作符对应,而且可以用于克隆。 BigInteger x = +y; 和 BigInteger x = y; 语句是完全不同的,前者为克隆一份新的对象,后者只是直接引用同一对象。之所以不实现 IClone 接口,是因为该接口的 Clone 方法的返回类型是 object,使用起来很不方便:var x = y.Clone() as BigInteger; 需要一次强制类型转换。如果 FCL 中有泛型版本的 IClone<T> 接口就好了。

 92     public static BigInteger operator ++(BigInteger x)
 93     {
 94       return x + 1;
 95     }
 96 
 97     public static BigInteger operator --(BigInteger x)
 98     {
 99       return x - 1;
100     }
101 

在 C# 语言中重载 ++ 和 -- 运算符比在 C++ 语言中容易,C# 编译器会处理好前缀和后缀的情况。

102     public static BigInteger operator +(BigInteger x, BigInteger y)
103     {
104       if (x == null || y == nullreturn null;
105       if (x.AbsCompareTo(y) < 0) Swap(ref x, ref y);
106       BigInteger z = new BigInteger(x);
107       if (x.sign * y.sign == -1) z.AbsSubtract(y);
108       else z.AbsAdd(y);
109       return z;
110     }
111 
112     public static BigInteger operator -(BigInteger x, BigInteger y)
113     {
114       if (x == null || y == nullreturn null;
115       return x + (-y);
116     }
117 
118     public static BigInteger operator *(BigInteger x, BigInteger y)
119     {
120       if (x == null || y == nullreturn null;
121       BigInteger z = 0;
122       z.sign = x.sign * y.sign;
123       z.data = new List<int>(new int[x.data.Count + y.data.Count]);
124       for (int i = x.data.Count - 1; i >= 0; i--)
125         for (int j = y.data.Count - 1; j >= 0; j--)
126         {
127           long n = Math.BigMul(x.data[i], y.data[j]);
128           z.data[i + j] += (int)(n % Base);
129           z.CarryUp(i + j);
130           z.data[i + j + 1+= (int)(n / Base);
131           z.CarryUp(i + j + 1);
132         }
133       z.Shrink();
134       return z;
135     }
136 
137     public static BigInteger operator /(BigInteger dividend, BigInteger divisor)
138     {
139       BigInteger remainder;
140       return DivRem(dividend, divisor, out remainder);
141     }
142 
143     public static BigInteger operator %(BigInteger dividend, BigInteger divisor)
144     {
145       BigInteger remainder;
146       DivRem(dividend, divisor, out remainder);
147       return remainder;
148     }
149 

在 C# 中重载加减乘除等运算符,C# 编译器也会自动添上相应的赋值运算符(+=, -=, ...)。

150     public static BigInteger DivRem(BigInteger dividend, BigInteger divisor, out BigInteger remainder)
151     {
152       remainder = null;
153       if (dividend == null || divisor == nullreturn null;
154       if (divisor.sign == 0throw new Exception("division by zero");
155       if (dividend.AbsCompareTo(divisor) < 0)
156       {
157         remainder = new BigInteger(dividend);
158         return 0;
159       }
160       remainder = 0;
161       BigInteger quotient = 0;
162       quotient.sign = dividend.sign * divisor.sign;
163       quotient.data = new List<int>(new int[dividend.data.Count]);
164       divisor = Abs(divisor); // NOT: divisor.sign = Math.Abs(divisor.sign);
165       for (int i = dividend.data.Count - 1; i >= 0; i--)
166       {
167         remainder = remainder * Base + dividend.data[i];
168         int iQuotient = remainder.BinarySearch(divisor, -1);
169         quotient.data[i] = iQuotient;
170         remainder -= divisor * iQuotient;
171       }
172       quotient.Shrink();
173       if (remainder.sign != 0) remainder.sign = dividend.sign;
174       return quotient;
175     }
176 
177     public static BigInteger Sqrt(BigInteger x)
178     {
179       if (x == null || x.sign < 0return null;
180       BigInteger root = 0;
181       root.sign = 1;
182       root.data = new List<int>(new int[x.data.Count / 2 + 1]);
183       for (int i = root.data.Count - 1; i >= 0; i--) root.data[i] = x.BinarySearch(root, i);
184       root.Shrink();
185       return root;
186     }
187 
188     int BinarySearch(BigInteger divisor, int i)
189     {
190       int low = 0, high = Base - 1, mid = 0, cmp = 0;
191       while (low <= high)
192       {
193         mid = (low + high) / 2;
194         cmp = CompareTo(divisor, mid, i);
195         if (cmp > 0) low = mid + 1;
196         else if (cmp < 0) high = mid - 1;
197         else return mid;
198       }
199       return (cmp < 0? (mid - 1) : mid;
200     }
201 
202     int CompareTo(BigInteger divisor, int mid, int i)
203     {
204       if (i >= 0) divisor.data[i] = mid;
205       return AbsCompareTo(divisor * ((i >= 0? divisor : mid));
206     }
207 
208     void AbsAdd(BigInteger x)
209     {
210       for (int i = 0; i < x.data.Count; i++)
211       {
212         data[i] += x.data[i];
213         CarryUp(i);
214       }
215     }
216 
217     void AbsSubtract(BigInteger x)
218     {
219       for (int i = 0; i < x.data.Count; i++)
220       {
221         data[i] -= x.data[i];
222         CarryDown(i);
223       }
224       Shrink();
225     }
226 
227     void CarryUp(int n)
228     {
229       for (; data[n] >= Base; n++)
230       {
231         if (n == data.Count - 1) data.Add(data[n] / Base);
232         else data[n + 1+= data[n] / Base;
233         data[n] %= Base;
234       }
235     }
236 
237     void CarryDown(int n)
238     {
239       for (; data[n] < 0; n++)
240       {
241         data[n + 1]--;
242         data[n] += Base;
243       }
244       Shrink();
245     }
246 
247     void Shrink()
248     {
249       for (int i = data.Count - 1; i >= 0 && data[i] == 0; i--) data.RemoveAt(i);
250       if (data.Count == 0) sign = 0;
251     }
252 

除法、取模和求平方根都是通过二分搜索来进行的,每次搜索时都有进行一次乘法和比较运算,以决定下次要搜索的区间。

253     public static bool operator ==(BigInteger x, BigInteger y)
254     {
255       if (object.ReferenceEquals(x, null)) return object.ReferenceEquals(y, null);
256       return x.Equals(y);
257     }
258 
259     public static bool operator !=(BigInteger x, BigInteger y)
260     {
261       if (object.ReferenceEquals(x, null)) return !object.ReferenceEquals(y, null);
262       return !x.Equals(y);
263     }
264 
265     public static bool operator <(BigInteger x, BigInteger y)
266     {
267       if (object.ReferenceEquals(x, null)) return !object.ReferenceEquals(y, null);
268       return x.CompareTo(y) < 0;
269     }
270 
271     public static bool operator >(BigInteger x, BigInteger y)
272     {
273       if (object.ReferenceEquals(x, null)) return false;
274       return x.CompareTo(y) > 0;
275     }
276 
277     public static bool operator <=(BigInteger x, BigInteger y)
278     {
279       if (object.ReferenceEquals(x, null)) return true;
280       return x.CompareTo(y) <= 0;
281     }
282 
283     public static bool operator >=(BigInteger x, BigInteger y)
284     {
285       if (object.ReferenceEquals(x, null)) return object.ReferenceEquals(y, null);
286       return x.CompareTo(y) >= 0;
287     }
288 
289     public override string ToString()
290     {
291       StringBuilder sb = new StringBuilder();
292       if (sign < 0) sb.Append('-');
293       sb.Append((data.Count == 0? 0 : data[data.Count - 1]);
294       for (int i = data.Count - 2; i >= 0; i--) sb.Append(data[i].ToString("D" + Len));
295       return sb.ToString();
296     }
297 
298     public override int GetHashCode()
299     {
300       int hash = sign;
301       foreach (int n in data) hash ^= n;
302       return hash;
303     }
304 
305     public override bool Equals(object other)
306     {
307       if (other == null || GetType() != other.GetType()) return false;
308       return Equals(other as BigInteger);
309     }
310 
311     public bool Equals(BigInteger other)
312     {
313       return CompareTo(other) == 0;
314     }
315 
316     public int CompareTo(BigInteger other)
317     {
318       if (object.ReferenceEquals(other, null)) return 1;
319       if (sign < other.sign) return -1;
320       if (sign > other.sign) return 1;
321       if (sign == 0return 0;
322       return sign * AbsCompareTo(other);
323     }
324 
325     int AbsCompareTo(BigInteger other)
326     {
327       if (data.Count < other.data.Count) return -1;
328       if (data.Count > other.data.Count) return 1;
329       for (int i = data.Count - 1; i >= 0; i--)
330         if (data[i] != other.data[i])
331           return (data[i] < other.data[i]) ? -1 : 1;
332       return 0;
333     }
334   }
335 }

最后就是各种逻辑运算符以及从基类继承的 ToString、GetHashCode、Equals 方法,以及用于实现接口的 Equals 和 CompareTo 方法。


下面,让我们来比较这四种 BigInteger 的效率吧。测试程序如下所示:

TestMain.cs:

 1 using System;
 2 
 3 namespace Skyiv
 4 {
 5   class TestMain
 6   {
 7     static void Main()
 8     {
 9       string s1 = new TestFcl35BigInteger().Run(20);
10       string s2 = new TestSkyivBigInteger().Run(20);
11       string s3 = new TestChewBigInteger().Run(20);
12       string s4 = new TestJavaBigInteger().Run(14);
13       Console.WriteLine(s1 == s2 && s1 == s3);
14     }
15   }
16 }

TestBase.cs:

 1 using System;
 2 using System.Diagnostics;
 3 
 4 namespace Skyiv
 5 {
 6   class TestBase
 7   {
 8     public string Run(int n)
 9     {
10       var stopWatch = Stopwatch.StartNew();
11       var s = RunTest(n);
12       stopWatch.Stop();
13       Console.WriteLine("{0,-25} {1,11:F7} {2,7:N0}",
14         GetType(), stopWatch.Elapsed.TotalSeconds, s.Length);
15       return s;
16     }
17 
18     protected virtual string RunTest(int n)
19     {
20       throw new InvalidOperationException("TestBase.RunTest()");
21     }
22   }
23 }

Test?????BigInteger.cs (以 TestSkyivBigInteger 为例,其它的大同小异):

 1 using System;
 2 using BigInteger = Skyiv.Numeric.BigInteger;
 3 
 4 namespace Skyiv
 5 {
 6   sealed class TestSkyivBigInteger : TestBase
 7   {
 8     protected sealed override string RunTest(int n)
 9     {
10       BigInteger a = 2, b;
11       for (int i = 2; i <= n; i++, a *= b) b = a + i;
12       return a.ToString();
13     }
14   }
15 }

这个测试程序主要是测试 BigInteger 的加法和乘法运算,其运行结果如下:

Skyiv.TestFcl35BigInteger  25.8100647 256,142
Skyiv.TestSkyivBigInteger  34.2497761 256,142
Skyiv.TestChewBigInteger  396.8586529 256,142
Skyiv.TestJavaBigInteger   61.9083458   4,003
True

上表中,第二栏是测试程序的运行时间,单位是秒。第三栏是最终运算结果的十进制数字的个数。前三个 BigInteger 都使用了同一测试用例,运算结果达 256,142 个十进制数字,并且经过比较这三个 BigInteger 的运算结果是相同的。最后的 java.math.BigInteger 由于运算速度极慢(估计是微软 J# 的问题,我想在真正的 java 语言中应该是比较快的),测试用例的规模更小,运算结果只有 4,003 个十进制数字。从测试结果可以看出:

  • Fcl35.Numeric.BigInteger 最快,推荐使用。
  • Skyiv.Numeric.BigInteger 也很快,并且源程序代码短小精悍,也推荐使用。
  • ChewKeongTAN.BigInteger 很慢,并且不能动态增长,不推荐使用。
  • java.math.BigInteger 极慢,并且需要 J# run-time library,也不推荐使用。

Skyiv.Numeric.BigInteger 和 ChewKeongTAN.BigInteger 的乘法都是使用简单的乘法,时间复杂度是 O(N2)。Fcl35.Numeric.BigInteger 的乘法是使用 Karatsuba 算法,其时间复杂度约为 O(N1.585),如下所示:

 1 public static BigInteger Multiply(BigInteger x, BigInteger y)
 2 {
 3   int xLength = x.Length;
 4   int yLength = y.Length;
 5   if ((((xLength + yLength) >= UpperBoundForSchoolBookMultiplicationDigits /* 0x40 */)
 6     && (xLength >= ForceSchoolBookMultiplicationThresholdDigits /* 8 */))
 7     && (yLength >= ForceSchoolBookMultiplicationThresholdDigits /* 8 */))
 8   {
 9     return MultiplyKaratsuba(x, y);
10   }
11   return MultiplySchoolBook(x, y);
12 }
13 
14 private static BigInteger MultiplyKaratsuba(BigInteger x, BigInteger y)
15 {
16   int length = Math.Max(x.Length, y.Length) / 2;
17   if (length <= 2 * ForceSchoolBookMultiplicationThresholdDigits /* 0x10 */
18     || x.Length < 2 * ForceSchoolBookMultiplicationThresholdDigits /* 0x10 */
19     || y.Length < 2 * ForceSchoolBookMultiplicationThresholdDigits /* 0x10 */)
20     return MultiplySchoolBook(x, y);
21   int shift = BitsPerDigit /* 0x20 */ * length;
22   BigInteger i1 = RightShift(x, shift);
23   BigInteger i2 = x.RestrictTo(length);
24   BigInteger i3 = RightShift(y, shift);
25   BigInteger i4 = y.RestrictTo(length);
26   BigInteger i5 = Multiply(i1, i3);
27   BigInteger i6 = Multiply(i2, i4);
28   BigInteger i7 = Multiply(i1 + i2, i3 + i4) - (i5 + i6);
29   return (i6 + LeftShift(i7 + LeftShift(i5, shift), shift));
30 }

其实如果使用快速傅里叶变换(FFT)来进行乘法运算,时间复杂度可降低到 O(N logN loglogN)。而且除法可以用除数的倒数乘以被除数来计算,倒数值用牛顿迭代法:

Ui+1 = Ui (2 - VUi)

来计算,这导致 U二次收敛于 1/V。实际上,许多超级计算机和 RISC 机都是使用这种迭代法来实现除法的。用牛顿迭代法计算平方根和除法类似。若:

Ui+1 = Ui (3 - VUi2) / 2

则 U二次收敛于 1/√V,最后用 V 除就得到√V

我可能会使用快速傅里叶变换来改写 Skyiv.Numeric.BigInteger 类。微软声称在下一版本的 .NET Framework 中会重新发布 System.Numeric.BigInteger,希望也使用快速傅里叶变换。:)

 

本文中的所有源代码可以到 http://bitbucket.org/ben.skyiv/biginteger/downloads 页面下载。
也可以使用 hg clone http://bitbucket.org/ben.skyiv/biginteger/ 命令下载。
关于 hg ,请参阅 Mercurial 备忘录

 


  1. 2008-07-25 更新: 请参阅:再谈 BigInteger - 使用快速傅里叶变换
  2. 2010-04-22 更新: 对 Skyiv.Numeric.BigInteger 的 GetHashCode、Parse 和 ToString 方法进行了优化。请参阅:再谈 BigInteger - 优化
posted on 2008-07-13 21:25 银河 阅读(...) 评论(...) 编辑 收藏