阶乘算法大全

发表于   |   分类于   |     |   阅读次数

简述

百度到的阶乘算法大全，先收藏一下，说不定有用。

原文：http://blog.chinaunix.net/uid-20788636-id-1841373.html

模板源代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
//阶乘各算法的C++类实现
#include <cmath>
#include <iomanip>
#include <cstring>
#include <iostream>
using namespace std;

class Factorial
{
    static const int MAXN = 5001; // 最大阶乘数，实际用不到这么大
    int *data[MAXN]; // 存放各个数的阶乘
    int *nonzero; // 从低位数起第一个非0数字
    int maxn; // 存放最大已经计算好的n的阶乘
    int SmallFact(int n); // n <= 12的递归程序
    void TransToStr(int n, int *s); // 将数n倒序存入数组中
    void Multply(int *A, int *B, int *C, int totallen);  // 执行两个高精度数的乘法

public:
    Factorial();
    ~Factorial();
    void Calculate(int n); // 调用计算阶乘
    int FirstNonZero(int n); // 返回阶乘末尾第一个非0数字
    int CountZeros(int n); // 返回阶乘末尾有多少个0
    int SecondNum(int n); // 返回阶乘左边的第二个数字
    bool CanDivide(int m, int n); // 判断数值 m 是否可以整除 n!
    void Output(int n) const;
};

int Factorial::SmallFact(int n)
{
    if (n == 1 || n == 0)
        return 1;
    return SmallFact(n - 1) * n;
}

void Factorial::TransToStr(int n, int *tmp)
{
    int i = 1;
    while (n)
    {
        tmp[i++] = n % 10;
        n /= 10;
    }
    tmp[0] = i - 1;
}

void Factorial::Multply(int *A, int *B, int *C, int totallen)
{
    int i, j, len;
    memset(C, 0, totallen * sizeof(int));
    for (i = 1; i <= A[0]; i++)
    {
        for (j = 1; j <= B[0]; j++)
        {
            C[i + j - 1] += A[i] * B[j]; // 当前i+j-1位对应项 + A[i] * B[j]
            C[i + j] += C[i + j - 1] / 10; // 它的后一位 + 它的商（进位）
            C[i + j - 1] %= 10; // 它再取余即可
        }
    }
    len = A[0] + B[0];
    while (len > 1 && C[len] == 0)
    {
        len--;    // 获得它的实际长度
    }
    C[0] = len;
}

Factorial::Factorial()   // 构造函数，先把<=12的阶乘计算好
{
    maxn = 12;
    data[0] = new int [2];
    data[0][0] = 1;
    data[0][1] = 1;
    int i, j = 1;
    for (i = 1; i <= 12; i++)
    {
        data[i] = new int [12];
        j = j * i;
        TransToStr(j, data[i]);
    }
    nonzero = new int [10 * MAXN];
    nonzero[0] = 1;
    nonzero[1] = 1; // nonzero[0]存储已经计算到的n!末尾非0数
}

Factorial::~Factorial()
{
    for (int i = 0; i <= maxn; i++)
        delete []data[i];
    delete []nonzero;
}

void Factorial::Calculate(int n)
{
    if (n > MAXN)
        return;
    if (n <= maxn)
        return;    // <= maxn的，已经在计算好的数组中了
    int i, j, len;
    int tmp[12];
    for (i = maxn + 1; i <= n; i++)
    {
        TransToStr(i, tmp);
        len = data[i - 1][0] + tmp[0] + 1;
        data[i] = new int [len + 1];
        Multply(data[i - 1], tmp, data[i], len + 1);
    }
    maxn = n;
}

int Factorial::FirstNonZero(int n)
{
    if (n >= 10 * MAXN)
    {
        cout << "Super Pig, your input is too large, cannot Calculate. Sorry! " << endl;
        return -1;
    }
    if (n <= nonzero[0])
        return nonzero[n];    //已经计算好了，直接返回

    int res[5][4] = {{0, 0, 0, 0}, {2, 6, 8, 4}, {4, 2, 6, 8}, {6, 8, 4, 2}, {8, 4, 2, 6}};
    int i, five, t;
    for (i = nonzero[0] + 1; i <= n; i++)
    {
        t = i;
        while (t % 10 == 0)
        {
            t /= 10;    // 先去掉 i 末尾的 0，这是不影响的
        }
        if (t % 2 == 0) // t是偶数直接乘再取模10即可
        {
            nonzero[i] = (nonzero[i - 1] * t) % 10;
        }
        else   // 否则转换成 res 数组来求
        {
            five = 0;
            while (t % 5 == 0)
            {
                if (five == 3)
                    five = 0;
                else
                    five++;
                t /= 5;
            }
            nonzero[i] = res[((nonzero[i - 1] * t) % 10) / 2][five];
            // (nonzero[i-1]*t)%10/2 正好序号为：1, 2, 3, 4 中的一个
        }
    }
    nonzero[0] = n;
    return nonzero[n];
}

/* 阶乘末尾有多少个0，实际上只与5的因子数量有关，即求 n/5+n/25+n/625+...... */
int Factorial::CountZeros(int n)
{
    if (n >= 2000000000)
    {
        cout << "Super Pig, your input is too large, cannot Calculate. Sorry! " << endl;
        return -1;
    }
    int cnt = 0;
    while (n)
    {
        n /= 5;
        cnt += n;
    }
    return cnt;
}

/* 输出N!左边第二位的数字：用实数乘，超过100就除以10，最后取个位即可 */
int Factorial::SecondNum(int n)
{
    if (n <= 3)
        return 0;
    int i;
    double x = 6;
    for (i = 4; i <= n; i++)
    {
        x *= i;
        while (x >= 100)
        {
            x /= 10;
        }
    }
    return (int(x)) % 10;
}

bool Factorial::CanDivide(int m, int n)
{
    if (m == 0)
        return false;
    if (n >= m)
        return true;
    int nn, i, j, nums1, nums2;
    bool ok = true;
    j = (int) sqrt(1.0 * m);
    for (i = 2; i <= j; i++)
    {
        if (m % i == 0)
        {
            nums1 = 0; // 除数m的素因子i的数量
            while (m % i == 0)
            {
                nums1++;
                m /= i;
            }
            nums2 = 0;
            nn = n;
            while (nn)   // 求 n 含有 i 因子的数量
            {
                nn /= i;
                nums2 += nn;
            }
            if (nums2 < nums1)   // 少于m中所含i的数量，则m肯定无法整除n!
            {
                ok = false;
                break;
            }
            j = (int)sqrt(1.0 * m); // 调整新的素因子前进范围
        }
    }
    if (!ok || m > n || m == 0)
        return false;
    else
        return true;
}

void Factorial::Output(int n) const
{
    if (n > MAXN)
    {
        cout << "Super Pig, your input is so large, cannot Calculate. Sorry! " << endl;
        return;
    }
    int i, len = 8;
    cout << setw(4) << n << "! = "; // 格式控制输出

    for (i = data[n][0]; i >= 1; i--)
    {
        cout << data[n][i];
        if (++len == 58)   // 实际每输出50个字符就换行
        {
            len = 8;
            cout << " ";
        }
    }
    if (len != 8)
        cout << endl;
}

int main()
{
    int n, m, i;
    Factorial f;
    while (cin >> n)
    {
        f.Calculate(n);
        f.Output(n);
        cout << "该阶乘末尾第一个非0数字是: " << f.FirstNonZero(n) << endl;
        cout << "该阶乘总共拥有数字0的个数：" << f.CountZeros(n) << endl;
        cout << "该阶乘的左边的第2位数字是：" << f.SecondNum(n) << endl;
        cin >> m;
        if (f.CanDivide(m, n))
            cout << m << " 可以整除 " << n << "! " << endl;
        else
            cout << m << " 不能整除 " << n << "! " << endl;
    }
    return 0;
}

高精度计算阶乘

这实际上是最没有技术含量的问题，但是又会经常用到，所以还是得编写，优化它的计算。

首先看小于等于12的阶乘计算（计算结果不会超出32位范围）：

1
2
3
4
5
int factorial(int n)
{
    if (n == 1 || n == 0) return 1;
    return factorial(n - 1) * n;
}

这个递归程序简单明了，非常直观，然而一旦n>12，则超过32位int型的范围出现错误结果，所以上面这个递归程序仅适合n<=12的阶乘计算。为了计算较大n的阶乘，需要将高精度乘法算法纳入到阶乘计算中来，高精度乘法过程可以如下简单的描述（其中A * B = C，A[0]，B[0]，C[0]分别存储长度）：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
for (i = 1; i <= A[0]; i++)
{
    for (j = 1; j <= B[0]; j++)
    {
        C[i + j - 1] += A[i] * B[j];        // 当前i+j-1位对应项 + A[i] * B[j]
        C[i + j] += C[i + j - 1] / 10;      // 它的后一位 + 它的商（进位）
        C[i + j - 1] %= 10;                 // 它再取余即可
    }
}
C[0] = A[0] + B[0];
while (C[0] > 1 && C[C[0]] == 0)
{
    C[0]--;   // 去头0，获得实际C的长度
}

有了这个高精度乘法之后，计算阶乘就可以简单的迭代进行：

1
2
3
4
5
for (i = 2; i <= n; i++)
{
    //将i转换成字符数组;
    //执行高精度乘法：将上一次结果乘上i
}

与数论有关

由于阶乘到后面越来越大，巧妙的利用数论求得一些有趣的数字（数值）等成为阶乘算法的设计点，下面给出几道相关的问题与分析：

计算阶乘末尾第一个非0数字：

这是一个比较经典的问题，比较复杂的算法是利用一个艰难的数学公式，可惜我不会，从网上的资料学习中，整理出下面这个简单易懂的算法：

观察n!，可以发现在乘的过程中，对于任意 n > 1，n!的末尾第一个非0数字都是偶数。我们只需保留最后一位非零数。当要乘的数中含有因数5时，我们可以把所有的因数5都当作8来乘。这是因为:

…x25=…10（舍）或…60，最后一位非零数为6。而恰好28=16，末位为6。

…x45=…70（舍）或…20，最后一位非零数为2。而恰好48=32，末位为2。

…x65=…30（舍）或…80，最后一位非零数为8。而恰好68=48，末位为8。

…x85=…90（舍）或…40，最后一位非零数为4。而恰好88=64，末位为4。

(对于n > 1时，最后一位不会出现 1, 7, 3, 9，而永远是2, 4, 6, 8的循环出现)

因此，在迭代作乘法时，主要就是计算因子5的数量，同时可见因子5的个数以4为循环节（即只需要取它的数量对4取模）。那么对于不同情况下的因子5的数量，可以通过

1
res[5][4] = {{0,0,0,0}, {2,6,8,4}, {4,2,6,8}, {6,8,4,2}, {8,4,2,6}};

来得到，使用nonzero[i]表示i的阶乘的最后一位，那么：

如果t是偶数，则直接乘：

1
nonzero[i] = (nonzero[i - 1] * t) % 10;

否则：

1
nonzero[i] = res[((nonzero[i - 1] * t) % 10) / 2][five];

其中t是除掉所有因子5的结果，five为因子5数量对4的模。相关题目：
http://acm.zju.edu.cn的第1222题。不过这一道题注意的是，它的输入n并非只在32位int数值范围内，而是有很大的长度，所以计算这道变态题目时，需要利用到高精度除法(n /= 5)和高精度加法(cnt += n)。

阶乘末尾有多少个0

分析发现，实际上形成末尾0，就是因子5的数量，而计算1~n之间包含一个因子i的个数的简单算法就是：

1
2
3
4
5
6
cnt = 0;
while (n)
{
    n /= i;
    cnt += n;
}

因此，直接将i换成5，就可以得到因子5的数量，也即n!末尾0的数量。相关题目：http://acm.zju.edu.cn的第2022题。

返回阶乘左边的第二个数字

简单算法：用实数乘，超过100就除以10，最后取个位即可。因为整数部分的个位就是阶乘结果左边的第二个数字。相关题目：http://acm.tongji.edu.cn的1016题。

判断数值m是否可以整除n!

算法：使用素因子判断法

A. 首先直接输出两种特殊情况：

m == 0 则0肯定不会整除n!；

n >= m 则m肯定可以整除n!;

B. 那么就只剩最后一种情况：m > n，我们从m的最小素因子取起，设素因子为i那么可以求得m的素因子i的个数nums1；再检查闭区间i ~ n之间的数，一共包含多少个素因子i，就可以简单的利用上面(2)中所介绍的数学公式进行计算得到nums2。如果nums2 < nums1，就表示1 ~ n中包含素因子的数量 < 除数m包含素因子i的数量，那么m必然不能整除n!，置ok = false。

C. 最后：如果!ok or m > n or m == 0则不能整除；否则可以整除

相关题目：http://acm.zju.edu.cn的第1850题。

数字N能否表示成若干个不相同的阶乘的和：

这里可以选择的阶乘为：0! ~ 9!，实际上这一题与数论无关，与搜索有关。相关题目：http://acm.zju.edu.cn的2358题。

分析，由于可供选择的阶乘数量较少，直接可以利用DFS搜索来做：

A. 首先将0 ~ 9的阶乘作一个表A[10]；再设置一个可以组成“和”的数组ans[N]。

B. 深度优先搜索方法：

1
2
3
4
5
6
7
8
9
10
search(n)
{
    for (i = n; i <= 9; i++)
    {
        sum += A[i];        //求和
        //如果sum在ans数组中不存在，则将sum插入到ans[]数组中
        search(n + 1);
        sum -= A[i];         //回溯
    }
}

C. 最后对于输入n，就在ans数组中查找是否存在n，如果存在，则表示n可以表示成不同的阶乘和，否则不行。