完全平方数的判定及整数平方根的快速求解

在上个星期的“有道难题”网络预赛中,某些Group的第二题涉及到了完全平方数的判定问题。相信大多数人的判定语句都差不多是这样写的:if(sqr((int)sqrtl(n))==n){...},把这条语句整理一下,大概是这样:

inline __int64 sqr( int n){ return  static_cast<__int64>(n)*n;} 

inline  int directSqrt( int n) 

     int root =  static_cast< int>(sqrtl(n)); //求n的开方,但是double类型的
     if (sqr(root+ 1) <= n) 
     { 
         return root +  1
     } 
     else  return root; 


inline  bool isSquareNumber( int n) 

   int root = directSqrt(n); 
   return sqr(root) == n; 
}

这个判定方法主要是利用标准库的sqrt函数来求整数平方根(integer square root),正确性不用怀疑,效率也足够。但是如果想要更快呢?通常来说,要判定任意一个数是否是完全平方数,求根的过程在所难免。但是,标准库的sqrt函数是求浮点平方根的,而我们只需要获取整数平方根即可,所以在这里肯定还有性能提升的空间。于是接下来就列举一些计算整数平方根的方法(后两个算法不是我写的),按照计算速度从慢到快依次排列,其中最后一个算法的效率大概是上述直接利用sqrtl的方法的四倍。

1. 线性搜索,直接从小到大依次取数试探,最原始的办法,代码如下:

inline  int linearSqrt( int n) 

     int i =  0
     while (sqr(i) <= n) ++i; 
     return i- 1
}


2. 二分法寻找根

inline  int bisectionSqrt( int n) 

     int l =  0, r = n/ 2 +  1
     while (r - l >  1
     { 
         int m = (l+r)/ 2
         __int64 sqrm = sqr(m); 
         if (sqrm == n) 
         { 
             return m; 
         } 
         else  if (sqrm > n) r = m; 
         else l = m; 
     } 
     
     if (sqr(r) == n) 
     { 
         return r; 
     } 
     else  return l; 
}

3. 牛顿迭代法,令方程 x^2 - n = 0,则求sqrt(n)也即是求这个方程的根,根据牛顿法的迭代公式有:x' = x - (x^2-n)/(2x) = (x^2+n)/(2x),选定合适的初始值迭代,很快就能收敛到根。另外,这个迭代其实也就是不动点迭代的松弛形式,由 x^2= n => x = n/x, 得不动点迭代x' = n/x,加上松弛系数1/2(算数平均),得x'= 1/2 * (x + n/x),即 x' = (x^2+n)/(2x),与前面牛顿法的迭代公式相同。算法的代码如下:

inline  int newtonSqrt( int n) 

     static  int table[ 16] = { 
         11<< 21<< 41<< 61<< 81<< 101<< 121<< 141<< 16
         1<< 181<< 201<< 221<< 241<< 261<< 281<< 30
     }; 
     unsigned i =  0, x0 =  0, x1 =  0
     while (table[i] < n && i!= 16) ++i; x0 =  1 << i; //寻找迭代的初识值
     while (x0 !=  0
     {     
         x1 = static_cast< int>((sqr(x0) + n)/ ( 2*x0)); 
         if (x1 >= x0) 
         { 
             break
         } 
         x0 = x1; 
     } 
     return x0; 
}

4. 逐位确定法,从高到低依次确定根的各个二进制位。见下面这个代码,看起来似乎很飘逸,但实际上算法很简单明晰,具体的解释可以搜索James Ulery的论文:"Computing Integer Square Roots".

inline  int fastIntSqrt( int n) 

     int temp, nHat =  0, b =  0x8000, bshft =  15
     do 
     { 
         if (n >= (temp = (((nHat<< 1)+b)<<bshft--))) 
         { 
             nHat += b; 
             n -= temp; 
         } 
     }  while (b >>=  1); 
     return nHat; 
}

5.建表+牛顿迭代,这个算法号称是要用来替代java数学包中的(int)(java.lang.Math.sqrt(integer))函数的,速度奇快,首先建立了一张强大的表来确定初值,然后最多经过两次牛顿迭代后就得到精确解。

inline  int javaIntSqrt( int n) 

     static  int table[ 256] = { 
         0,     16,   22,   27,   32,   35,   39,   42,   45,   48,   50,   53,   55,   57
         59,    61,   64,   65,   67,   69,   71,   73,   75,   76,   78,   80,   81,   83
         84,    86,   87,   89,   90,   91,   93,   94,   96,   97,   98,   99101102
         103104106107108109110112113114115116117118
         119120121122123124125126128128129130131132
         133134135136137138139140141142143144144145
         146147148149150150151152153154155155156157
         158159160160161162163163164165166167167168
         169170170171172173173174175176176177178178
         179180181181182183183184185185186187187188
         189189190191192192193193194195195196197197
         198199199200201201202203203204204205206206
         207208208209209210211211212212213214214215
         215216217217218218219219220221221222222223
         224224225225226226227227228229229230230231
         231232232233234234235235236236237237238238
         239240240241241242242243243244244245245246
         246247247248248249249250250251251252252253
         253254254255 
     }; 

     int xn; 
     if (n >=  0x7FFEA810return  0xB504
     if (n >=  0x10000) { 
         if (n >=  0x1000000) { 
             if (n >=  0x10000000) { 
                 if (n >=  0x40000000) { 
                     xn = table[n>>  24] <<  8
                 }  else { 
                     xn = table[n>>  22] <<  7
                 } 
             }  else { 
                 if (n>=  0x4000000) { 
                     xn = table[n>>  20] <<  6
                 }  else { 
                     xn = table[n>>  18] <<  5
                 } 
             } 

             xn = (xn +  1 + (n/ xn)) >>  1
             xn = (xn +  1 + (n/ xn)) >>  1
             return ((xn * xn) > n) ? --xn : xn; 
         }  else { 
             if (n>=  0x100000) { 
                 if (n>=  0x400000) { 
                     xn = table[n>>  16] <<  4
                 }  else { 
                     xn = table[n>>  14] <<  3
                 } 
             }  else { 
                 if (n>=  0x40000) { 
                     xn = table[n>>  12] <<  2
                 }  else { 
                     xn = table[n>>  10] <<  1
                 } 
             } 

             xn = (xn +  1 + (n/ xn)) >>  1
             return ((xn * xn) > n) ? --xn : xn; 
         } 
     }  else { 
         if (n>=  0x100) { 
             if (n>=  0x1000) { 
                 if (n>=  0x4000) { 
                     xn = (table[n>>  8]) +  1
                 }  else { 
                     xn = (table[n>>  6] >>  1) +  1
                 } 
             }  else { 
                 if (n>=  0x400) { 
                     xn = (table[n>>  4] >>  2) +  1
                 }  else { 
                     xn = (table[n>>  2] >>  3) +  1
                 } 
             } 
             return ((xn * xn) > n) ? --xn : xn; 
         }  else { 
             if (n>=  0) { 
                 return table[n] >>  4
             } 
         } 
     } 
     return - 1
}

依次计算从1到1000000的每个数的整数平方根,下面分别是上面6个算法所用的时间:

0. 利用标准库sqrtl, directSqrt : 610 ms
1. 线性搜索, linearSqrt : 81256 ms
2. 二分法, bisectionSqrt : 2734 ms
3. 牛顿迭代法, newtonSqrt : 691 ms
4. 逐位确定法, fastIntSqrt : 301 ms
5. 建表+牛顿迭代, javaIntSqrt : 150 ms

可见牛顿迭代法与利用标准库sqrtl的方法速度相当,逐位确定法的速度大概是牛顿迭代法的两倍。当然,最快的还是建表+牛顿迭代,看来以其来替代(int)(java.lang.Math.sqrt(integer))也不是没有道理的。另,以上测试数据相对来说都比较小,可以想象,若是计算接近于2^31这样的大数的整数平方根,逐位确定法和建表+牛顿迭代应该更占据优势。

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值