fork(4) download
  1. #include <bits/stdc++.h>
  2. using namespace std;
  3. const int N = 1e8;
  4. bool comp[N];
  5. vector<int> primes;
  6.  
  7. void normal_sieve() {
  8. memset(comp, 0, sizeof comp);
  9.  
  10. primes.clear();
  11. bool done = false;
  12. for (int i = 2; i < N; i++) {
  13. if (!comp[i]) {
  14. primes.push_back(i);
  15.  
  16. if (i * i > N) {
  17. done = true;
  18. }
  19.  
  20. if (!done)
  21. for (int j = i * i; j < N; j += i)
  22. comp[j] = 1;
  23. }
  24. }
  25. }
  26.  
  27. void normal_optimized_sieve() {
  28. memset(comp, 0, (sizeof comp) / 2);
  29. primes.clear();
  30.  
  31. primes.push_back(2);
  32. int n = N / 2;
  33. int sr = round(sqrt(N));
  34. for (int i = 1; i < n; i++) {
  35. if (!comp[i]) {
  36. int p = 2 * i + 1;
  37. primes.push_back(p);
  38.  
  39. if (p > sr)
  40. continue;
  41.  
  42. for (int j = (p + 1) * i; j < n; j += p)
  43. comp[j] = 1;
  44. }
  45. }
  46. }
  47.  
  48. void linear_sieve() {
  49. memset(comp, 0, sizeof comp);
  50. primes.clear();
  51. for (int i = 2; i < N; i++) {
  52. if (!comp[i])
  53. primes.push_back(i);
  54. for (int j = 0, si = primes.size(); j < si && i * primes[j] < N; j++) {
  55. comp[primes[j] * i] = 1;
  56. if (i % primes[j] == 0)
  57. break;
  58. }
  59. }
  60. }
  61.  
  62. void linear_sieve_optimized() {
  63. memset(comp, 0, (sizeof comp) / 2);
  64. primes.clear();
  65. primes.push_back(2);
  66. for (int i = 1, p = 3; p < N; i++, p+=2) {
  67. if (!comp[i])
  68. primes.push_back(p);
  69. for (int j = 1, si = primes.size(); j < si && p * primes[j] < N; j++) {
  70. comp[(primes[j] * p) >> 1] = 1;
  71. if (p % primes[j] == 0)
  72. break;
  73. }
  74. }
  75. }
  76.  
  77. void segmented_sieve() {
  78. int n = N;
  79. const int S = 50000;
  80.  
  81. primes.clear();
  82. int nsqrt = round(sqrt(n));
  83. vector<char> is_prime(nsqrt + 1, true);
  84.  
  85. ::primes.push_back(2);
  86. vector<char> block(S);
  87.  
  88. vector<int> primes, start_indices;
  89. for (int i = 3; i <= nsqrt; i += 2) {
  90. if (is_prime[i]) {
  91. primes.push_back(i);
  92. start_indices.push_back(i * i);
  93. for (int j = i * i; j <= nsqrt; j += 2 * i)
  94. is_prime[j] = false;
  95. }
  96. }
  97.  
  98. for (int k = 0; k * S <= n; k++) {
  99. fill(block.begin(), block.end(), true);
  100. int start = k * S;
  101. for (auto i = 0u; i < primes.size(); i++) {
  102. int p = primes[i];
  103. int idx = start_indices[i];
  104. for (; idx < S; idx += 2 * p)
  105. block[idx] = false;
  106. start_indices[i] = idx - S;
  107. }
  108. if (k == 0)
  109. block[1] = false;
  110. for (int i = 1; i < S && start + i <= n; i += 2) {
  111. if (block[i])
  112. ::primes.push_back(start + i);
  113. }
  114. };
  115. }
  116.  
  117. void segmented_sieve_optimized() {
  118. int n = N;
  119. const int S = 30000;
  120.  
  121. primes.clear();
  122. int nsqrt = round(sqrt(n));
  123.  
  124. vector<char> is_prime(nsqrt + 1, true);
  125. vector<int> primes, start_indices;
  126. for (int i = 3; i <= nsqrt; i += 2) {
  127. if (is_prime[i]) {
  128. primes.push_back(i);
  129. start_indices.push_back((i * i - 1) / 2);
  130. for (int j = i * i; j <= nsqrt; j += 2 * i)
  131. is_prime[j] = false;
  132. }
  133. }
  134.  
  135. ::primes.push_back(2);
  136. vector<char> block(S);
  137. int high = (n - 1) / 2;
  138. for (int low = 0; low <= high; low += S) {
  139. fill(block.begin(), block.end(), true);
  140. for (auto i = 0u; i < primes.size(); i++) {
  141. int p = primes[i];
  142. int idx = start_indices[i];
  143. for (; idx < S; idx += p)
  144. block[idx] = false;
  145. start_indices[i] = idx - S;
  146. }
  147. if (low == 0)
  148. block[0] = false;
  149. for (int i = 0; i < S && low + i <= high; i++) {
  150. if (block[i])
  151. ::primes.push_back((low + i) * 2 + 1);
  152. }
  153. };
  154. }
  155.  
  156. void aitken() {
  157. int n = N;
  158. std::vector<bool> m(n + 1, false);
  159.  
  160. int root = std::sqrt(n) + 1;
  161. for (int i = 1; i <= root; i++) {
  162. for (int j = 1; j <= root; j++) {
  163. int a = 4 * i * i + j * j;
  164. int b = 3 * i * i + j * j;
  165. int c = 3 * i * i - j * j;
  166.  
  167. if (a <= n && (a % 12 == 1 || a % 12 == 5))
  168. m[a].flip();
  169. if (b <= n && b % 12 == 7)
  170. m[b].flip();
  171. if (i > j && c <= n && c % 12 == 11)
  172. m[c].flip();
  173. }
  174. }
  175.  
  176. for (int i = 5; i * i <= n; i++) {
  177. if (m[i]) {
  178. for (int j = 1; j * i * i <= n; j++)
  179. m[j * i * i] = false;
  180. }
  181. }
  182.  
  183. primes.clear();
  184. primes.push_back(2);
  185. primes.push_back(3);
  186. for (int i = 5; i < n; i++) {
  187. if (m[i])
  188. primes.push_back(i);
  189. }
  190. }
  191.  
  192. template <typename T>
  193. void time_compare(string name, T func, vector<int> comparison) {
  194. primes.clear();
  195. auto t0 = clock();
  196. (*func)();
  197. auto t1 = clock();
  198. if (!comparison.empty())
  199. assert(comparison == primes);
  200. cout << setw(20) << name << ": " << (1.0 * (t1 - t0)) / CLOCKS_PER_SEC << " seconds"
  201. << endl;
  202. }
  203.  
  204. int main() {
  205. ::primes.reserve(6'000'00);
  206.  
  207. time_compare("normal", normal_sieve, {});
  208. auto comparison = primes;
  209. time_compare("normal_optimized", normal_optimized_sieve, comparison);
  210. time_compare("linear", linear_sieve, comparison);
  211. time_compare("linear_optimized", linear_sieve_optimized, comparison);
  212. time_compare("segmented", segmented_sieve, comparison);
  213. time_compare("segmented_optimized", segmented_sieve_optimized, comparison);
  214. time_compare("aitken", aitken, comparison);
  215. }
  216.  
Success #stdin #stdout 3.6s 147648KB
stdin
Standard input is empty
stdout
              normal: 1.09295 seconds
    normal_optimized: 0.468311 seconds
              linear: 0.546933 seconds
    linear_optimized: 0.321325 seconds
           segmented: 0.188974 seconds
 segmented_optimized: 0.133987 seconds
              aitken: 0.800168 seconds