Miller_Rabin大素数测试 与 Pollard_rho整数分解模版

1. #include <cstdio>  
2. #include <cstring>  
3. #include <cstdlib>  
4. #include <algorithm>  
5. using namespace std;  
6. typedef __int64 LL;  
7. const int Times = 20;  
8. LL factor[100], l;  
9. LL gcd(LL a, LL b)  
10. {  
11.     return b ? gcd(b, a%b):a;  
12. }  
13. LL add_mod(LL a, LL b, LL n)  
14. {  
15.     LL ans = 0;  
16.     while(b)  
17.     {  
18.         if(b&1)  
19.             ans = (ans + a)%n;  
20.         b >>= 1;  
21.         a = (a<<1)%n;  
22.     }  
23.     return ans;  
24. }  
25. LL pow_mod(LL a, LL m, LL n)  
26. {  
27.     LL ans = 1;  
28.     while(m)  
29.     {  
30.         if(m&1)  
31.             ans = add_mod(ans, a, n);  
32.         m >>= 1;  
33.         a = add_mod(a, a, n);  
34.     }  
35.     return ans;  
36. }  
37. bool Witness(LL a, LL n)  
38. {  
39.     int j = 0;  
40.     LL m = n-1;  
41.     while(!(m&1))  
42.     {  
43.         j++;  
44.         m >>= 1;  
45.     }  
46.     LL x = pow_mod(a, m, n);  
47.     if(x == 1 || x == n-1)  
48.         return true;  
49.     while(j--)  
50.     {  
51.         x = add_mod(x, x, n);  
52.         if(x == n-1)  
53.             return true;  
54.     }  
55.     return false;  
56. }  
57. bool Miller_Rabin(LL n)  
58. {  
59.     if(n < 2)  
60.         return false;  
61.     if(n == 2)  
62.         return true;  
63.     if(!(n&1))  
64.         return false;  
65.     for(int i = 0; i < Times; i++)  
66.     {  
67.         LL a = rand()%(n-1)+1;  
68.         if(!Witness(a, n))  
69.             return false;  
70.     }  
71.     return true;  
72. }  
73. LL Pollard_rho(LL n, LL c)  
74. {  
75.     LL i = 1, x = rand()%(n-1)+1, y = x, k = 2, d;  
76.     //srand(time(NULL));  
77.     while(true)  
78.     {  
79.         i++;  
80.         x = (add_mod(x,x,n)+c)%n;  
81.         d = gcd(y-x,n);  
82.         if(d > 1 && d < n)  
83.             return d;  
84.         if(y == x)  
85.             return n;  
86.         if(i == k)  
87.         {  
88.             y = x;  
89.             k <<= 1;  
90.         }  
91.     }  
92. }  
93. void get_fact(LL n, LL k)  
94. {  
95.     if(n == 1)  
96.         return;  
97.     if(Miller_Rabin(n))  
98.     {  
99.         factor[l++] = n;  
100.         return;  
101.     }  
102.     LL p = n;  
103.     while(p >= n)  
104.     {  
105.         p = Pollard_rho(p, k--);  
106.     }  
107.     get_fact(p, k);  
108.     get_fact(n/p, k);  
109. }