（转载）求解阶乘乘积的最右非零值

转载自：http://hi.baidu.com/shilyx/blog/item/5bb7733e6313ec3e70cf6cd9.html，感谢作者分享！

I.问题的提出
整数n的阶乘写作n!，它是从l到n的所有整数的乘积。阶乘增长的速度很快：13!在大多数计算机上不能用32位的整数来存放。70!已经超出多数浮点类型的范围。你的任务是找出n!最右边的非零位。例如，5!＝1*2*3*4*5＝120，所以5!的最右非零位为2，同样，7!＝1*2*3*4*5*6*7＝5040，所以7!的最右非零位为4。50!=30414093201713378043612608166064768844377641568960512000000000000，所以最右非零位为2。
写一个程序，计算N(1<=N<=50,000,000)阶乘的最右边的非零位的值。注意:10,000,000！有2499999个零。
输入：一个正整数n。
输出：n!的最右非零位。
函数接口：int factorial(int n);

为表述方便，可以以“尾数”的概念指代最右非零值这一说法。

II.问题分析
对于这么一个问题，恐怕没有人真的把阶乘值算出来，稍大的数据的阶乘值就非常大了。既然是去求最后一位，那么我们可以只运算最后一位，其余数位无须过问。我们可以从数字1开始，依次做乘法运算，遇到零则舍去，然后继续考虑上一个数位。在这种思想下，可以轻松得到一个C语言的算法：
int factorial(int n)
{
int result = 1;
for(int i = 1; i <= n; i++)
{
result *= i;
while(result % 10 == 0)
{
result /= 10;
}
result %= 10;
}
return result;
}
代入数据50检验，得输出为4，这个答案是错误的。跟踪程序的运行并和正确结果做对比，可以发现在运算到数字15的时候首次出现了错误。我们来分析一下。
14！=87178291200
15!= 15*14!=1307674368000
按照我们刚才的分析，在运算到14!之后，我们在程序中保存了一个2，这是14!的最后一位，我们用2乘以15得30，这时再舍弃零，得结果3，但是真实的结果是8。很明显的，12*15和2*15结果是不样的。因此这个算法不正确，之所以会出现错误，是因为最后一位的运算中产生了一个整十，因此要向前一位找到新的非零数，而我们没有保留这一位的全部信息，因此将得到错误的结果。

III.整十进位产生的根源
对于一般的进位(不是整十)我们是不需要特别考虑的，这些进位发生以后，仍然保存着最后一位非零，不会影响我们的运算结果。那么整十进位是怎么产生的呢？
大于十的数计算阶乘的时候，中间就直接夹杂着10这一因子，它可以导致整十进位。所有的因子5同因子2乘起来，也得到10，数字10本身也是5和2的乘积。因此将整个乘式因子分解后，我们可以断定：成对出现的2和5是产生10的根源。单独的2或5和一切其它数字的乘法，都不能产生整十。我们可以先把所有的5*2从乘式中剔除，然后再使用上面提到的第一种算法计算结果就可以了。
然而在乘式中，5的个数和2的个数未必是相等的，实际上除了在1!和0!中两者都不存在因子2和因子5外，在其它乘式中两者数量根本不能相等，后面我们会证明这一点。按照直觉，我们会猜测2更多一点，所以在剔除了等量的因子5和因子2之后，因子5就不存在了，但是因子2还剩下一部分，我们只需要处理剩下的数。但是究竟是谁剩下呢？

IV.5和2，谁更关键？
由于5>2因此因子5的个数更少，因子5的个数决定了整十的个数，因此也决定了乘式结尾零的个数。下面我们设计一算法来计算因子5的个数：

/*计算对Number!做因子分解后Factor因子的个数并返回*/
int GetFactorNumber(unsigned Number, unsigned Factor)
{
int result = 0;
if(Factor < 2 || Number < 1)
return -1;
while(Number >= Factor)
result += Number /= Factor;
return result;
}

V.计算最右非零位
按照上面的思路，我们现在要着手计算阶乘结果的尾数。对于一个给定的数字M，我们要求得M!中的因子5的个数N，然后从1到M依次处理数据，并在这些数据中除去所有的因子5和N个因子2，然后使用余下的数据按照最开始时提出的算法计算结果。代码及相关分析如下：

/*计算Number!的最右非零位并输出*/
int GetLastDigit(unsigned Number)
{
int result = 1;/*这是连乘法的基础值，也是最后要返回的结果*/
int i;/*循环变量*/
int tmp;/*循环中用到的临时值*/
int Count = 0;/*记录找到的因子2的个数，直到和因子5的个数相等*/
int TotalFactors = GetFactorNumber(Number, 5);/*因子5的总数*/
for(i = Number; i >= 1; i--)/*从Number一直向下处理*/
{
tmp = i;
/*尽可能多的得到因子2，直到tmp是个奇数或者不需要更多的因子2*/
while(tmp % 2 == 0 && Count < TotalFactors)
{
tmp /= 2;
Count++;
}
while(tmp % 5 == 0)/*除去所有的因子5*/
tmp /= 5;
tmp %= 10;/*取得个位数*/
result *= tmp;/*做阶乘*/
result %= 10;/*只保留个位数*/
}
return result;/*返回运算结果*/
}

代入数据检验，输出结果完全吻合。
其实不仅仅是最右一位数可以求得，最右任意位数都可以按照上面的方法计算，只是要保留的数据位数不同而已。

VI.更深的探索
上面找到了计算阶乘积最后几位数的一般方法，这个方法对每一个数据都进行了处理，因此，速度上受到了一点限制。考虑阶乘的计算方法细节，我们可以发现一些有用的东西。
(1). 阶乘结果的尾数只能是2,4,8之中的一个(除0!和1!外)
为什么会这样呢？我们可以这样来思考，对于7以下的数字，我们可以非常容易的验证这一点，就不赘述了；对于7及以上数字的阶乘，里面必然包含一个非素因子6，我们可以把阶乘结果表示6A，A是其余数字的乘积，设6A的最后一位数字是k。那么我们考虑6K*6这个数值，且不论这个数值有什么意义，我们先推断一下这个数值的特征。它的尾数是什么呢？6A*6=6*6A，而6*6=36的尾数依然是6，因此6A*6的尾数等同于6A的尾数，也是k。这也就是说，6*k的尾数也是k本身。
在数字1-9中，满足这一条的只有三个数字:2,4,8。
(2).不同数字对阶乘尾数的影响有差异
考虑数字1和6，它们对尾数的变化是没有影响的，尾数与其相乘之后保持不变。
考虑数字9和4，它们对尾数的变化有影响，但是影响是周期性的。尾数同两个9相乘，相当于与1相乘，因为9*9=91；尾数同两个4相乘，相当于与6相乘，因为4*4=16。
考虑数字2，3，7，8，它们对尾数的变化有影响，但是影响也是周期性的。尾数同四个2或者四个8相乘，相当于与6相乘；尾数同四个3或者四个7相乘，相当于与1相乘；
数字0和5无须考虑，所有的5的倍数都已经被完全分解。
根据上面的条件分析，可以得到一个新的思路，这个思路是基于为尾数分布的分析的。比如对尾数9来说，如果阶乘中含有奇数个数字9那么就相当于一个数字9，偶数个数字9相当于一个数字1。
(3).统计尾数的可行性
设法统计到尾数出现的次数，我们的工作会更加轻松。我们现在要确定的是从数字1到N，其间尾数为k的数有多少(k=1..9)。这种分布是以10为周期的，在N和与N接近的一个整十之间，也可能存在以k为尾数的数。那么我们可以根据下面是算法来统计这些数的数目：

/* 统计1-Number之间的数字中以Digit为尾数的个数*/
int GetLastDigitTimes(unsigned Number, int Digit)
{
if(Number % 10 >= Digit)
return Number / 10 + 1;
return Number / 10;
}

VII.更优的算法
基于以上的分析，我们可以设计一个新的算法，此算法不要求将每个数据处理一次，因此，时间上会有一定的节省。我们现在可以一开始就统计数据吗？还不可以。因为除掉因子5和一部分因子2的工作还没有做，这部分工作我没有找到更好的办法，因此还要使用上面的算法，直到除去了足够多的因子2，然后再使用新算法。对于M的新算法的流程如下：
(1).计算因子5的个数；
(2).从M向1依次处理数据，尽可能多的除去因子2，直到除去的足够多；
(3).从下一个数据开始，使用新算法。计算余下数字中的2,3,4,7,8,9的个数，并分别处理；
(4).分析统计结果，得出最终结果。
算法及其分析如下：

/* 计算Number的Times次方*/
int Power(int Number, int Times)
{
int result = 1;
while(Times--)
result *= Number;
return result;
}

/* 计算Number!的最右非零位并输出方法2*/
int GetLastDigitB(unsigned Number)
{
int result = 1;/*这是连乘法的基础值，也是最后要返回的结果*/
int i;/*循环变量*/
int tmp;/*循环中用到的临时值*/
int Count = 0;/*记录找到的因子2的个数，直到和因子5的个数相等*/
int TotalFactors = GetFactorNumber(Number, 5);/*因子5的总数*/
int T2, T3, T4, T7, T8, T9;/*指示以2,3,4,7,8,9为尾数的个数*/
/*下面依旧是排除因子2的过程*/
for(i = Number; i >= 1; i--)
{
tmp = i;
while(tmp % 2 == 0 && Count < TotalFactors)
{
tmp /= 2;
Count++;
}
while(tmp % 5 == 0)
tmp /= 5;
tmp %= 10;
result *= tmp;
result %= 10;
if(Count == TotalFactors && i % 5 == 1)
break;
}
i--;

/*因子2排除完毕，开始使用新方法统计各尾数出现次数并对自身的周期取模*/
T2 = GetLastDigitTimes(i, 2) % 4;
T3 = GetLastDigitTimes(i, 3) % 4;
T4 = GetLastDigitTimes(i, 4) % 2;
T7 = GetLastDigitTimes(i, 7) % 4;
T8 = GetLastDigitTimes(i, 8) % 4;
T9 = GetLastDigitTimes(i, 9) % 2;

/*处理余下数据中的所有因子5*/
for(; i >= 1; i -= 5)
{
tmp = i;
while(tmp % 5 == 0)
tmp /= 5;
tmp %= 10;
result *= tmp;
result %= 10;
}

/*处理得到的不同尾数的出现次数*/
result *= Power(2, T2);
result *= Power(3, T3);
result *= Power(4, T4);
result %= 10;
result *= Power(7, T7);
result %= 10;
result *= Power(8, T8);
result %= 10;
result *= Power(9, T9);
result %= 10;
return result;
}

这个新的算法在很大程度上减小了循环体的长度，如果能找到更好的处理因子2的办法，那么算法还可以更优。

VIII.两种算法的比较
主要是从运行时间是比较两种算法，在Number值较小的情况下，不好判断哪一个速度更快，我们计算一亿的阶乘的尾数，并使用Windows提供的GetTickCount函数来标志运行时间，在我的机器上前一个算法耗费时间是10484，新的算法则是5016，由此可见，新的算法在大量数据处理上可以节约大约一半的时间。