素数

1. 普通筛........

2. 埃氏筛

　　复杂度听说是O(nloglogn) 听说证明挺鬼畜

 1 #include <bits/stdc++.h>
 2 using namespace std;
 3 
 4 const int MAXN = 11111111;
 5 bool isPrime[MAXN];
 6 
 7 int main(){
 8     //freopen("test.txt", "w", stdout);
 9     memset(isPrime, true, sizeof(isPrime));
10     int x = 0, k = ceil(sqrt(MAXN));
11     for(int i = 2; i < k; ++i){
12         if(isPrime[i])
13             for(int j = i; i * j < MAXN; ++j)
14                 isPrime[i*j] = false;
15     }
16     for(int i = 2; i < MAXN; ++i) 
17         if(isPrime[i])
18             printf("%d ", i);
19     return 0;
20 }

3. 欧拉筛

　　还是跟埃氏筛一样用算术基本定理

　　但是埃氏筛对于某些数会多次筛除

　　欧拉筛使合数只被其最小质因数筛除时间复杂度为 O(n)

 1 #include <bits/stdc++.h>
 2 using namespace std;
 3 
 4 const int MAXN = 11111111;
 5 int prime[MAXN];
 6 bool isPrime[MAXN];
 7 
 8 int main(){
 9     //freopen("book.txt", "w", stdout);
10     memset(isPrime, true, sizeof(isPrime));
11     int x = 0;
12     for(int i = 2; i < MAXN; ++i){
13         if(isPrime[i]) prime[x++] = i;
14         for(int j = 0; j < x; ++j){
15             if(i * prime[j] >= MAXN) break;
16             isPrime[i*prime[j]] = false;
17             if(i % prime[j] == 0) break;
18         }
19     }
20     for(int i = 0; i < x; ++i) printf("%d ", prime[i]);
21     return 0;
22 }

精髓是 if(i % prime[j] == 0) break;

当某个数能被某质数整除时，后面的这个数的倍数的最小质因子一定是当前的质数而不是其他的质数，所以退出循环

例如：当 i == 6 时筛去(6 * 2 = )12后，6被2整除，说明6*3、6*5的最小质因数是2，而不是3和5

4.Miller_rabin强伪素数测试

Int64以内Rabin-Miller强伪素数测试和Pollard 因数分解的算法实现

 1 typedef long long ll;
 2 ll mul_mod(ll a, ll b, ll m){
 3     ll ans = 0, exp = a;
 4     while(exp >= m) exp -= m;
 5     while(b){
 6         if(b & 1){
 7             ans += exp;
 8             while(ans >= m) ans -= m;
 9         }
10         exp += exp;
11         while(exp >= m) exp -= m;
12         b >>= 1;
13     }
14     return ans;
15 }
16 ll pow_mod(ll x, ll n, ll m){
17     ll ans = 1, exp = x;
18     while(exp >= m) exp -= m;
19     while(n){
20         if(n & 1) ans = mul_mod(ans, exp, m);
21         exp = mul_mod(exp, exp, m);
22         n >>= 1;
23     }
24     return ans;
25 }
26 bool miller_rabin(ll n){
27     if(n < 2) return false;
28     if(n == 2) return true;
29     if(!(n & 1)) return false;
30     ll p = n - 1;
31     int k = 0;
32     while(!(p & 1)){
33         ++k;
34         p >>= 1;
35     }
36     for(int i = 0; i < 20; ++i){
37         ll a = rand() % (n - 1) + 1;
38         ll x = pow_mod(a, p, n);
39         if(x == 1) continue;
40         bool flag = false;
41         for(int j = 0; j < k; ++j){
42             if(x == n - 1){
43                 flag = true;
44                 break;
45             }
46             x = mul_mod(x, x, n);
47         }
48         if(flag) continue;
49         return false;
50     }
51     return true;
52 }

附：POJ 2429 代码

  1 #include <map>
  2 #include <vector>
  3 #include <iostream>
  4 #include <algorithm>
  5 using namespace std;
  6 
  7 typedef long long ll;
  8 
  9 ll mul_mod(ll a, ll b, ll mod){
 10     ll ans = 0, exp = a;
 11     while(exp >= mod) exp -= mod;
 12     while(b){
 13         if(b & 1){
 14             ans += exp;
 15             if(ans >= mod) ans -= mod;
 16         }
 17         exp <<= 1;
 18         if(exp >= mod) exp -= mod;
 19         b >>= 1;
 20     }
 21     return ans;
 22 }
 23 ll pow_mod(ll x, ll n, ll mod){
 24     ll ans = 1, exp = x;
 25     while(exp >= mod) exp -= mod;
 26     while(n){
 27         if(n & 1) ans = mul_mod(ans, exp, mod);
 28         exp = mul_mod(exp, exp, mod);
 29         n >>= 1;
 30     }
 31     return ans;
 32 }
 33 
 34 const int MP = 233333;
 35 vector<int> primes;
 36 vector<bool> isprime;
 37 void init_primes(){
 38     isprime = vector<bool>(MP, true);
 39     isprime[0] = isprime[1] = false;
 40     for(int i = 2; i < MP; ++i){
 41         if(isprime[i]) primes.push_back(i);
 42         for(int j = 0; j < primes.size() && i * primes[j] < MP; ++j){
 43             isprime[i * primes[j]] = false;
 44             if(i % primes[j] == 0) break;
 45         }
 46     }
 47 }
 48 bool miller_rabin(ll n){
 49     if(n < 2) return false;
 50     if(n == 2) return true;
 51     if(!(n & 1)) return false;
 52     ll p = n - 1;
 53     int k = 0;
 54     while(!(p & 1)){
 55         ++k;
 56         p >>= 1;
 57     }
 58     for(int i = 0; i < 20; ++i){
 59         ll a = rand() % (n - 1) + 1;
 60         ll x = pow_mod(a, p, n);
 61         if(x == 1) continue;
 62         bool flag = false;
 63         for(int j = 0; j < k; ++j){
 64             if(x == n - 1){
 65                 flag = true;
 66                 break;
 67             }
 68             x = mul_mod(x, x, n);
 69         }
 70         if(flag) continue;
 71         return false;
 72     }
 73     return true;
 74 }
 75 bool is_prime(ll n){
 76     if(n < MP) return isprime[n];
 77     else return miller_rabin(n);
 78 }
 79 
 80 ll pollard_rho(ll n, int a){
 81     ll x = 2, y = 2, d = 1;
 82     while(d == 1){
 83         x = mul_mod(x, x, n) + a;
 84         y = mul_mod(y, y, n) + a;
 85         y = mul_mod(y, y, n) + a;
 86         if(x == y) break;
 87         d = __gcd((x >= y ? x - y : y - x), n);
 88     }
 89     if(d == 1 || d == n) return pollard_rho(n, a + 1);
 90     return d;
 91 }
 92 void factorize(ll n, map<ll, int>& factors){
 93     if(is_prime(n)) ++factors[n];
 94     else{
 95         for(int i = 0; i < primes.size(); ++i)
 96             while(!(n % primes[i])){
 97                 ++factors[primes[i]];
 98                 n /= primes[i];
 99             }
100         if(n != 1){
101             if(is_prime(n)) ++factors[n];
102             else{
103                 ll d = pollard_rho(n, 1);
104                 factorize(d, factors);
105                 factorize(n / d, factors);
106             }
107         }
108     }
109 }
110 
111 int main(){
112     ios::sync_with_stdio(false);
113     cin.tie(0);
114     ll a, b;
115     init_primes();
116     while(cin >> a >> b){
117         ll c = b / a;
118         map<ll, int> factors;
119         factorize(c, factors);
120         vector<ll> p_factors;
121         for(map<ll, int>::iterator it = factors.begin(); it != factors.end(); ++it){
122             ll tmp = 1;
123             while(it->second--) tmp *= it->first;
124             p_factors.push_back(tmp);
125         }
126         ll x = 1, y = c, sum = 1e18;
127         for(int i = 0; i < (1 << p_factors.size()); ++i){
128             ll tx = 1;
129             for(int j = 0; j < p_factors.size(); ++j)
130                 if(i & (1 << j))
131                     tx *= p_factors[j];
132             ll ty = c / tx;
133             if(tx < ty && tx + ty < sum){
134                 x = tx;
135                 y = ty;
136                 sum = x + y;
137             }
138         }
139         cout << x * a << " " << y * a << endl;
140     }
141     return 0;
142 }