AtCoder NoviStepsを埋めてみる(7) 再帰全探索 1Q以降の続きです。今回はエラトステネスの篩です。
Contents
エラトステネスの篩とは？
エラトステネスの篩は、N 以下の正の素数を O(N loglog N) の計算量で列挙するアルゴリズムです。古代ギリシアの科学者、エラトステネスが考案したといわれています。
1 ～ n までの数字を書く
1 は素数ではないので消す
残った数のなかで最小の 2 は素数。それ以外の 2 の倍数は素数ではないので消す。
残った数のなかで次に最小の 3 は素数。それ以外の 3 の倍数は素数ではないので消す。
残った数のなかで次に最小の 5 は素数。それ以外の 5 の倍数は素数ではないので消す。
これを繰り返すと最終的に素数のみが残る。
篩い落とし操作は n の平方根 までやれば充分です。これだと n の平方根 より大きな素数による合成数を消すことができませんが、これは n よりも大きな値となるので気にする必要はありません。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 |
// エラトステネスの篩で max 以下の整数が素数かどうかを調べる static bool[] Eratosthenes(int max) { bool[] is_primes = new bool[max + 10]; Array.Fill(is_primes, true); is_primes[0] = false; is_primes[1] = false; int ceiling = (int)(Math.Sqrt(max)) + 1; for (int cur = 2; cur <= ceiling; cur++) { if (!is_primes[cur]) continue; for (long v = 1L * cur * cur; v <= max; v += cur) is_primes[v] = false; } return is_primes; } |
あとこんなメソッドも定義しておきます。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 |
// エラトステネスの篩の結果から max 以下のすべての素数のリストを返す static List<int> GetPrimes(int max) { List<int> res = new List<int>(); var is_primes = Eratosthenes(max); for (int i = 0; i <= max; i++) { if (is_primes[i]) res.Add(i); } return res; } |
B26 – Output Prime Numbers
問題の概要
N 以下の素数を値の小さい方から全て出力せよ。
これは基本問題です。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 |
class Program { // bool[] Eratosthenes(int max) // List<int> GetPrimes(int max) 上記定義のとおり static void Main() { int X = int.Parse(Console.ReadLine()); List<int> primes = GetPrimes(X); foreach (int v in primes) { Console.WriteLine(v); } } } |
D – 2017-like Number
問題の概要
N も (N + 1) ÷ 2 も素数であるという条件を満たす奇数 N を 2017に似た数 と定義する。
Q 個のクエリで l[i], r[i] が与えられる。
l[i] 以上 r[i] 以下の 2017に似た数 となる奇数 の個数を求めよ。
ある奇数が 2017に似た数かどうかはエラトステネスの篩で得た配列を見ればわかるのですが、l[i] 以上 r[i] 以下の 2017に似た数 となる奇数 の個数をクエリのたびに数えていては時間がかかります。
そこで前処理として 10^5以下の整数について 2017 に似た数かどうかを調べるとともに sums[i] に i 以下の 2017 に似た数の個数を格納しています。クエリがきたら sums[r] – sums[l – 1] を計算することですぐに解を返せるようにしています。
|
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 |
class Program { // bool[] Eratosthenes(int max) // List<int> GetPrimes(int max) 上記定義のとおり static void Main() { int max = 100000; var is_primes = Eratosthenes(max + 10); int cnt = 0; int[] sums = new int[max + 10]; for (int i = 1; i < max; i++) { if (i % 2 == 1 && is_primes[i] && is_primes[(i + 1) / 2]) cnt++; sums[i] = cnt; } (int, int) ReadInt2() { string[] vs = Console.ReadLine().Split(); return (int.Parse(vs[0]), int.Parse(vs[1])); } int Q = int.Parse(Console.ReadLine()); for (int i = 0; i < Q; i++) { (int l, int r) = ReadInt2(); Console.WriteLine(sums[r] - sums[l - 1]); } } } |
D – 250-like Number
問題の概要
k が素数 p < q であるふたつの素数を使って k = p × q^3 と表されるとき、k を「 250 に似た数」と定義する。 N 以下の 250 に似た数 は全部でいくつあるだろうか？ N が 10^18 と大きいのですが、q が取りうる範囲も 10^6 なので全探索可能です。 q を全探索します。N を q^3 で割ってみてその整数部分と q - 1 を比較します。両者のうち小さいほうが p の最大値となります。 エラトステネスの篩で 10^6 以下の素数を全取得するとともに i 以下の素数の個数をすぐに求めることができるように counts[i] に i 以下の素数の個数を格納しておきます。
|
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 |
class Program { // bool[] Eratosthenes(int max) 上記定義のとおり static void Main() { int max = 1000000; // i 以下の素数の個数を求める bool[] is_primes = Eratosthenes(max + 10); List<int> primes = new List<int>(); int[] counts = new int[max + 10]; for (int i = 0; i <= max; i++) { if (is_primes[i]) primes.Add(i); counts[i] = primes.Count; } long N = long.Parse(Console.ReadLine()); long ans = 0; foreach (long q in primes) { if (N / (q * q * q) == 0) break; long a = N / (q * q * q); long b = q - 1; long min = Math.Min(a, b); ans += counts[min]; } Console.WriteLine(ans); } } |
D – Coprime 2
問題の概要
長さ N の正整数列 A が与えられる。
すべての要素と互いに素である 1 以上 M 以下の整数をすべて求めてください。
A のすべての要素がもたない素因数と、これらを掛け合わせて生成される合成数が出力すべき解です。
A のそれぞれの要素を素因数分解しようとすると時間がかかりそうですが、工夫すると短時間で終わります。
1 から M までの整数で取得できた素因数で割ることができるかどうかの処理もひとつひとつ確認していては時間がかかるのでエラトステネスの篩のような処理をすることで時間短縮しています。
まず素因数分解をするためのクラスを定義します。
|
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 |
class PrimeFactorization { List<long> Primes; HashSet<long> PrimeSet; public PrimeFactorization(int prime_max) { Primes = GetPrimes(prime_max + 10); PrimeSet = new HashSet<long>(Primes); } bool[] Eratosthenes(int max) { bool[] is_primes = new bool[max + 10]; Array.Fill(is_primes, true); is_primes[0] = false; is_primes[1] = false; int ceiling = (int)(Math.Sqrt(max)) + 1; for (int cur = 2; cur <= ceiling; cur++) { if (!is_primes[cur]) continue; for (long v = 1L * cur * cur; v <= max; v += cur) is_primes[v] = false; } return is_primes; } List<long> GetPrimes(int max) { List<long> res = new List<long>(); var is_primes = Eratosthenes(max); for (int i = 0; i <= max; i++) { if (is_primes[i]) res.Add(i); } return res; } public Dictionary<long, int> Solve(long v) { Dictionary<long, int> dic = new Dictionary<long, int>(); if (PrimeSet.Contains(v)) { dic.Add(v, 1); return dic; } foreach (long p in Primes) { while (v % p == 0) { v /= p; if (!dic.ContainsKey(p)) dic.Add(p, 0); dic[p]++; } if (PrimeSet.Contains(v)) { dic.Add(v, 1); v = 1; break; } if (v == 1) break; if (p * p > v) { dic.Add(v, 1); v = 1; break; } } // 割り切れずに残ったものは Primes.Max() よりも大きい素数である if (v != 1) dic.Add(v, 1); return dic; } // ついでに素数かどうか判定するメソッドも追加する public bool IsPrime(long v) { if(v <= 1) return false; if (PrimeSet.Contains(v)) return true; foreach (long p in Primes) { if (v % p == 0) return false; if (p * p > v) break; } return true; } } |
|
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 |
class Program { static void Main() { string[] vs = Console.ReadLine().Split(); int N = int.Parse(vs[0]); int M = int.Parse(vs[1]); int[] A = Console.ReadLine().Split().Select(_ => int.Parse(_)).ToArray(); int max = 100000 + 10; PrimeFactorization pf = new PrimeFactorization(max); HashSet<long> factors = new HashSet<long>(); // 素因数 foreach (int v in A) { var dic = pf.Solve(v); foreach (var pair in dic) factors.Add(pair.Key); } bool[] B = new bool[M + 1]; Array.Fill(B, true); foreach (long factor in factors) { for (long i = factor; i <= M; i += factor) B[i] = false; } List<int> ans = new List<int>(); for (int i = 1; i <= M; i++) { if (B[i]) ans.Add(i); } Console.WriteLine(ans.Count); foreach (var v in ans) { Console.WriteLine(v); } } } |
E – Coprime
問題の概要
長さ N の数列 A が与えられます。
① すべての 0 ≦ i, j ≦ N – 1 について、GCD(A[i] ,A[j]) = 1 が成り立つか？
② ① が成り立たないとき、GCD(A[0], A[i], …, A[N – 1]) = 1 は成り立つか？
① が成り立つときは pairwise coprime、② が成り立つときは setwise coprime、どちらでもない場合は not coprime と出力せよ。
① である場合はそれぞれの要素が共通する素因数を持たないことを意味しています。なので、それを確認すればよいです。② かどうかは実際に最大公約数を計算すればよいです。
最大公約数を求める処理はユークリッドの互除法を用いています。
|
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 |
// PrimeFactorization クラスは上記定義のとおり class Program { static void Main() { int N = int.Parse(Console.ReadLine()); int[] A = Console.ReadLine().Split().Select(_ => int.Parse(_)).ToArray(); int max = 1000000; PrimeFactorization pf = new PrimeFactorization(max); // それぞれの要素が共通する素因数を持たないことを確認する bool IsPairwiseCoprime(int[] A) { HashSet<long> factors = new HashSet<long>(); foreach (int v in A) { var dic = pf.Solve(v); foreach (var pair in dic) { if (!factors.Contains(pair.Key)) factors.Add(pair.Key); else return false; } } return true; } // すべての要素の最大公約数を求める bool IsSetwiseCoprime(int[] A) { int GCD(int a, int b) => b > a ? GCD(b, a) : (b == 0 ? a : GCD(b, a % b)); int gcd = A[0]; for (int i = 1; i < A.Length; i++) { gcd = GCD(gcd, A[i]); } return gcd == 1; } if (IsPairwiseCoprime(A)) Console.WriteLine("pairwise coprime"); else if (IsSetwiseCoprime(A)) Console.WriteLine("setwise coprime"); else Console.WriteLine("not coprime"); } } |
D – 9 Divisors
問題の概要
N 以下の正整数のうち、正の約数をちょうど 9 個持つものの個数を求めよ。
素因数分解したときに A^a × B^b × C^c × … となる自然数の約数の個数は (a + 1) × (b + 1) × (c + 1) × … です。
では (a + 1) × (b + 1) × (c + 1) × … が 9 になる場合とはどのような場合でしょうか？ a = 8 のときと、a = 2, b = 2 のとき以外ありません。
a = 8 のときは比較的簡単です。エラトステネスの篩で素数を取得してそれらの 8 乗が N を超えるまで数え上げます。a = 2, b = 2 の場合は実は同じ素数の組み合わせをダブルカウントしないように 2 つの素数 p, q を p < q と考えます。素数の 2 乗のリストを生成して、すべての p について N / (p * p) 以下で (p * p) + 1 以上の個数を二分探索法で数え上げます。p が大きくなると N / (p * p) が p を下回ったり、N / (p * p) 以下で (p * p) + 1 以上の要素の個数が 0 になります。その場合はそれ以上調べても見つからないので探索を終了します。
二分探索の処理は ライブラリ ac-library-csharp の StlFunction.BinarySearch メソッドを使用しています。
|
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 |
using AtCoder; class Program { // bool[] Eratosthenes(int max) // List<int> GetPrimes(int max) 上記定義のとおり static void Main() { int max = 2000000 + 10; List<int> primes = GetPrimes(max); long N = long.Parse(Console.ReadLine()); long ans = 0; // 素数の 8 乗パターン foreach (long p in primes) { long v = 1; for (int i = 0; i < 8; i++) v *= p; if (v > N) break; ans ++; } // 素数の 2 乗 × 素数の 2 乗 パターン List<long> prime2s = primes.Select(_ => 1L * _ * _).ToList(); for (int i = 0; i < prime2s.Count; i++) { long n = N; long p = prime2s[i]; n /= p; if (n < p) break; int ok = prime2s.Count; int ng = -1; int idx = StlFunction.BinarySearch(ok, ng, mid => prime2s[mid] > n); if (idx - i - 1 >= 0) ans += idx - i - 1; else break; } Console.WriteLine(ans); } } |
D – AABCC
問題の概要
N 以下の正整数のうち、a < b < c なる 素数 a, b, c を用いて a^2 × b × c^2 と表せるものはいくつあるか求めよ。
N が 10^12 と大きいですが、a < b < c という条件があるため、a を大きくしたらそれにともなって c も大きくするしかなく、意外と探索範囲は狭いです。
|
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 |
class Program { // bool[] Eratosthenes(int max) // List<int> GetPrimes(int max) 上記定義のとおり static void Main() { int max = 1000000 + 10; List<int> primes = GetPrimes(max); long N = long.Parse(Console.ReadLine()); HashSet<Int128> set = new HashSet<Int128>(); for (int x = 0; x < primes.Count; x++) { // a = primes[x] と決めた段階で a < b < c の条件より // a * a * b * c * c の値はこれよりも大きくなってしまう if (Int128.One * primes[x] * primes[x] * primes[x] * primes[x] * primes[x] > N) break; for (int y = x + 1; y < primes.Count; y++) { // b = primes[y] と決めた段階で a < b < c の条件より // a * a * b * c * c の値はこれよりも大きくなる if (Int128.One * primes[x] * primes[x] * primes[y] * primes[y] * primes[y] > N) break; for (int z = y + 1; z < primes.Count; z++) { Int128 v = Int128.One * primes[x] * primes[x] * primes[y] * primes[z] * primes[z]; if (v > N) break; set.Add(v); } } } Console.WriteLine(set.Count); } } |
097 – Primes in an Interval
問題の概要
L 以上 R 以下の素数の個数を求めよ。
これはエラトステネスの区間篩です。L と R が 10^12 に近い大きな数ですが、これらの数が合成数であれば 10^6 の素因数をもつため、エラトステネスの篩で 10^6 までの素数を求めて、この倍数をふるい落とすことで素数の個数を得ることは可能です。
まずエラトステネスの区間篩のためのクラスを定義します。
|
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 |
class PrimeInterval { bool[] IsPrimes; long Min; long Max; public PrimeInterval(long min, long max) { Min = min; Max = max; IsPrimes = new bool[Max - Min + 1]; Array.Fill(IsPrimes, true); List<long> primes = GetPrimes((int)Math.Sqrt(Max) + 10); foreach (long p in primes) { long start = Min % p == 0 ? Min / p * p : Min / p * p + p; for (long i = start; i <= Max; i += p) { if (i != p) IsPrimes[i - Min] = false; // 素数の倍数であれば合成数（ただし1倍は除く） } } // 上の処理だけでは 1 がふるい落とせていないので注意 if (Min == 0) IsPrimes[1] = false; if (Min == 1) IsPrimes[0] = false; } bool[] Eratosthenes(int max) { bool[] is_primes = new bool[max + 10]; Array.Fill(is_primes, true); is_primes[0] = false; is_primes[1] = false; int ceiling = (int)(Math.Sqrt(max)) + 1; for (int cur = 2; cur <= ceiling; cur++) { if (!is_primes[cur]) continue; for (long v = 1L * cur * cur; v <= max; v += cur) is_primes[v] = false; } return is_primes; } List<long> GetPrimes(int max) { List<long> res = new List<long>(); var is_primes = Eratosthenes(max); for (int i = 0; i <= max; i++) { if (is_primes[i]) res.Add(i); } return res; } public bool IsPrime(long v) { return IsPrimes[v - Min]; } } |
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 |
class Program { static void Main() { string[] vs = Console.ReadLine().Split(); long L = long.Parse(vs[0]); long R = long.Parse(vs[1]); PrimeInterval primeInterval = new PrimeInterval(L, R); int ans = 0; for (long i = L; i <= R; i++) { if (primeInterval.IsPrime(i)) ans++; } Console.WriteLine(ans); } } |
