数論・素数に関するメソッドを C# で一通り実装してまとめました。
通常の数値計算プログラミング用にも、また競技プログラミング用のライブラリとしても使えます。
機能は次の通りです。
Euclid 互除法
- 最大公約数 (GCD)
- 最小公倍数 (LCM)
ax + by = 1の解 (Euclid 互除法の拡張)
素数
下の O(x) は計算量のオーダーを表します。
- 素因数分解
- O(√n)
- 約数の列挙
- O(√n)
- 素数判定
- O(√n)
- n 以下の素数の列挙
- およそ O(n)
- m 以上 M 以下の素数の列挙
- およそ O(√M + (M – m))
ソースコードは以下の通りです。
ここではシンプルな実装とし、例外処理も全体的に省略しています。
さらに無駄な演算処理を除くように改良することで、もう少し速くすることはできるはずです。
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| using System; | |
| namespace AlgorithmLib.Numerics | |
| { | |
| public static class Euclidean | |
| { | |
| public static int Gcd(int a, int b) { for (int r; (r = a % b) > 0; a = b, b = r) ; return b; } | |
| public static int Lcm(int a, int b) => a / Gcd(a, b) * b; | |
| public static long Gcd(long a, long b) { for (long r; (r = a % b) > 0; a = b, b = r) ; return b; } | |
| public static long Lcm(long a, long b) => a / Gcd(a, b) * b; | |
| // ax + by = 1 の解 | |
| // 前提: a と b は互いに素。 | |
| // ax + by = GCD(a, b) の解を求める場合、予め GCD(a, b) で割ってからこの関数を利用します。 | |
| public static (long x, long y) ExtendedEuclid(long a, long b) | |
| { | |
| if (b == 0) throw new ArgumentOutOfRangeException(nameof(b)); | |
| if (b == 1) return (1, 1 - a); | |
| var q = Math.DivRem(a, b, out var r); | |
| var (x, y) = ExtendedEuclid(b, r); | |
| return (y, x - q * y); | |
| } | |
| } | |
| } |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| using System; | |
| using System.Collections.Generic; | |
| namespace AlgorithmLib.Numerics | |
| { | |
| public static class Primes | |
| { | |
| /// <summary> | |
| /// 指定された自然数を素因数分解します。O(√n) | |
| /// </summary> | |
| /// <param name="n">自然数。</param> | |
| /// <returns>素因数のコレクション。n = 1 の場合は空の配列。</returns> | |
| public static long[] Factorize(long n) | |
| { | |
| var r = new List<long>(); | |
| for (long x = 2; x * x <= n && n > 1; ++x) | |
| while (n % x == 0) | |
| { | |
| r.Add(x); | |
| n /= x; | |
| } | |
| // √n を超える素因数はたかだか 1 個であり、その次数は 1。 | |
| if (n > 1) r.Add(n); | |
| return r.ToArray(); | |
| } | |
| /// <summary> | |
| /// 指定された自然数の約数をすべて求めます。O(√n) | |
| /// </summary> | |
| /// <param name="n">自然数。</param> | |
| /// <returns>約数のコレクション。</returns> | |
| public static long[] Divisors(long n) | |
| { | |
| var r = new List<long>(); | |
| for (long x = 1; x * x <= n; ++x) | |
| if (n % x == 0) r.Add(x); | |
| for (int i = r.Count - 1; i >= 0; --i) | |
| { | |
| var v = n / r[i]; | |
| if (v > r[i]) r.Add(v); | |
| } | |
| return r.ToArray(); | |
| } | |
| /// <summary> | |
| /// 指定された自然数が素数かどうかを判定します。O(√n) | |
| /// </summary> | |
| /// <param name="n">自然数。</param> | |
| /// <returns>素数である場合に限り true。</returns> | |
| public static bool IsPrime(long n) | |
| { | |
| for (long x = 2; x * x <= n; ++x) | |
| if (n % x == 0) return false; | |
| return n > 1; | |
| } | |
| /// <summary> | |
| /// 指定された自然数以下の素数をすべて求めます。およそ O(n) | |
| /// </summary> | |
| /// <param name="n">自然数。</param> | |
| /// <returns>素数のコレクション。</returns> | |
| public static long[] GetPrimes(long n) | |
| { | |
| var b = new bool[n + 1]; | |
| for (long p = 2; p * p <= n; ++p) | |
| if (!b[p]) | |
| for (long x = p * p; x <= n; x += p) | |
| b[x] = true; | |
| var r = new List<long>(); | |
| for (long x = 2; x <= n; ++x) if (!b[x]) r.Add(x); | |
| return r.ToArray(); | |
| } | |
| /// <summary> | |
| /// 指定された範囲内の素数をすべて求めます。およそ O(√M + (M - m)) | |
| /// </summary> | |
| /// <param name="m">最小値。</param> | |
| /// <param name="M">最大値。</param> | |
| /// <returns>素数のコレクション。</returns> | |
| public static long[] GetPrimes(long m, long M) | |
| { | |
| var rM = (int)Math.Ceiling(Math.Sqrt(M)); | |
| var b1 = new bool[rM + 1]; | |
| var b2 = new bool[M - m + 1]; | |
| if (m == 1) b2[0] = true; | |
| for (long p = 2; p <= rM; ++p) | |
| if (!b1[p]) | |
| { | |
| for (var x = p * p; x <= rM; x += p) b1[x] = true; | |
| for (var x = Math.Max(p, (m + p - 1) / p) * p; x <= M; x += p) b2[x - m] = true; | |
| } | |
| var r = new List<long>(); | |
| for (int x = 0; x < b2.Length; ++x) if (!b2[x]) r.Add(m + x); | |
| return r.ToArray(); | |
| } | |
| } | |
| } |
使用例を次に示します。
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| using System; | |
| using AlgorithmLib.Numerics; | |
| namespace AlgorithmConsole | |
| { | |
| class Program | |
| { | |
| static void Main(string[] args) | |
| { | |
| // 素因数分解 | |
| Console.WriteLine(string.Join(" * ", Primes.Factorize(2020))); | |
| // 2 * 2 * 5 * 101 | |
| // 約数の列挙 | |
| Console.WriteLine(string.Join(" ", Primes.Divisors(2020))); | |
| // 1 2 4 5 10 20 101 202 404 505 1010 2020 | |
| // 素数判定 | |
| Console.WriteLine(Primes.IsPrime(1000000007)); | |
| // True | |
| // n 以下の素数の列挙 | |
| Console.WriteLine(string.Join(" ", Primes.GetPrimes(100))); | |
| // 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 | |
| // m 以上 M 以下の素数の列挙 | |
| Console.WriteLine(string.Join(" ", Primes.GetPrimes(1000000000, 1000000100))); | |
| // 1000000007 1000000009 1000000021 1000000033 1000000087 1000000093 1000000097 | |
| } | |
| } | |
| } |
前回: 二分探索のライブラリ化
作成したサンプル
バージョン情報
- C# 7.0
- .NET Standard 2.0
Do Design Space 