#1066804 (C#) No.8030 ミラー・ラビン素数判定法のテスト

提出ソース

結果

問題	No.8030 ミラー・ラビン素数判定法のテスト
ユーザー	crimsontea
提出日時	2025-04-12 02:09:53
言語	C# (.NET 8.0.404)
結果	AC
実行時間	1,127 ms / 9,973 ms
コード長	22,682 bytes
コンパイル時間	10,889 ms
コンパイル使用メモリ	170,140 KB
実行使用メモリ	201,952 KB
最終ジャッジ日時	2025-04-12 02:10:10
合計ジャッジ時間	15,721 ms
ジャッジサーバーID （参考情報）	judge2 / judge3

このコードへのチャレンジ
（要ログイン）

ファイルパターン	結果
other	AC * 10

権限があれば一括ダウンロードができます

コンパイルメッセージ

  復元対象のプロジェクトを決定しています...
  /home/judge/data/code/main.csproj を復元しました (122 ミリ秒)。
/home/judge/data/code/Main.cs(173,32): warning CS8632: '#nullable' 注釈コンテキスト内のコードでのみ、Null 許容参照型の注釈を使用する必要があります。 [/home/judge/data/code/main.csproj]
/home/judge/data/code/Main.cs(174,30): warning CS8632: '#nullable' 注釈コンテキスト内のコードでのみ、Null 許容参照型の注釈を使用する必要があります。 [/home/judge/data/code/main.csproj]
  main -> /home/judge/data/code/bin/Release/net8.0/main.dll
  main -> /home/judge/data/code/bin/Release/net8.0/publish/

ソースコード

raw source code

using System;
using System.Linq;
using System.Collections;
using System.IO;
using System.Collections.Generic;
using System.Text;
using System.Numerics;
using System.Runtime.Intrinsics.X86;
using System.Buffers;
using System.Diagnostics;
using System.Runtime.CompilerServices;
using CrimsonTeaLib;
using static CrimsonTeaLib.InputUtility;

class Program
{
    static void Main()
    {
        using Output output = new(false);

        Stopwatch sw = Stopwatch.StartNew();

        InputNewLine();
        var n = NextInt32;

        for (int i = 0; i < n; i++)
        {
            InputNewLine();
            var x = (ulong)NextInt64;

            var res = MillerRabinPrimalityTester.IsPrime(x) ? 1 : 0;
            Console.WriteLine($"{x} {res}");
        }

        Console.Error.WriteLine($"{sw.Elapsed.TotalMilliseconds} ms");
    }
}

public static class MillerRabinPrimalityTester
{
    private static readonly ulong[] Witnesses = { 2, 325, 9375, 28178, 450775, 9780504, 1795265022 };

    /// <summary>
    /// 指定された64ビット符号なし整数が素数かどうかを判定します。
    /// </summary>
    /// <param name="n">判定する整数。</param>
    /// <returns>nが素数の場合は true、それ以外の場合は false。</returns>
    /// <remarks>
    /// このメソッドはミラーラビン素数判定法を使用します。
    /// 計算量は O(k * (log n)^2) または O(k * (log n)^3) です (kはwitnessの数)。
    /// 64bitの範囲では、選択されたwitnessにより決定的な結果を返します。
    /// .NET 7.0 以降が必要です (UInt128 を使用するため)。
    /// </remarks>
    public static bool IsPrime(ulong n)
    {
        // 基本的なケースの処理
        if (n < 2) return false;
        if (n == 2 || n == 3) return true;
        if (n % 2 == 0) return false; // 2以外の偶数は素数ではない

        // 小さな奇数素数での試行除算 (効率化のため)
        // これらの素数で割り切れるなら、合成数であることが早くわかる
        if (n < Witnesses[^1]) // witness より小さい数は witness 自身と比較する必要がある
        {
            // witness配列に含まれるかチェック (小さい素数用)
            for (int i = 0; Witnesses[i] * Witnesses[i] <= n; ++i)
            {
                if (n == Witnesses[i]) return true;
                if (n % Witnesses[i] == 0) return false;
            }
            // 上記ループでカバーされない小さい素数のチェック
            // 例: 5, 7, 11, 13 ... (Witnessesに含まれないもの)
            // 必要であれば、より多くの小さな素数をここでチェックすることも可能
            // ただし、ミラーラビン本体でチェックされるため必須ではない
            if (n == 5 || n == 7 || n == 11 || n == 13 || n == 17 || n == 19 || n == 23 || n == 29 || n == 31 || n == 37) return true;
            if (n % 3 == 0 || n % 5 == 0 || n % 7 == 0 || n % 11 == 0 || n % 13 == 0 || n % 17 == 0 || n % 19 == 0 || n % 23 == 0 || n % 29 == 0 || n % 31 == 0 || n % 37 == 0) return false;
        }


        // ミラーラビンテストの準備
        // n - 1 = d * 2^s となる d と s を見つける
        ulong d = n - 1;
        int s = BitOperations.TrailingZeroCount(d); // 末尾の0の数を数える (dを2で割れる回数s)
        d >>= s; // d を 2^s で割る (dは奇数になる)

        // 各witnessでテストを実行
        foreach (ulong a in Witnesses)
        {
            // a % n == 0 の場合、n は a の倍数 (aが素数ならnは合成数、ただし n=a の場合を除く)
            // ここでは a >= n のチェックは不要 (witnessはnより小さい必要があるというより、a%nで考える)
            if (a % n == 0) continue; // 実質的に a >= n のケースを含む (n=aの場合を除く)

            // x = a^d mod n を計算
            ulong x = ModPow(a, d, n);

            // 条件1: x == 1
            // 条件2: x == n - 1 (ループの初回)
            if (x == 1 || x == n - 1) continue; // このwitnessでは素数の可能性がある

            bool probablePrime = false;
            // 条件2: x = a^(d*2^r) mod n == n - 1 (0 < r < s)
            for (int r = 1; r < s; r++)
            {
                x = ModMultiply(x, x, n); // x = x^2 mod n

                // 最適化: xが1になった場合
                // 前のステップで x != 1 かつ x != n-1 であり、現在のステップで x*x mod n = 1 となった。
                // これは n が合成数であることを示す (1の非自明な平方根が存在するため)。
                if (x == 1) return false;

                if (x == n - 1)
                {
                    probablePrime = true; // 条件2を満たした
                    break;
                }
            }

            // ループを抜けても probablePrime が false のままなら、
            // a^d mod n != 1 かつ、全ての r (0 <= r < s) について a^(d*2^r) mod n != n-1
            // これは n が合成数であることを示す (フェルマーテスト + 平方根テストに失敗)
            if (!probablePrime) return false;
        }

        // すべてのwitnessでテストを通過した -> 素数である (このwitnessセットでは確定的)
        return true;
    }

    /// <summary>
    /// (a * b) mod m を計算します。ulongのオーバーフローを防ぎます。
    /// .NET 7.0 以降が必要です。
    /// </summary>
    private static ulong ModMultiply(ulong a, ulong b, ulong m)
    {
        // UInt128 を使用して 128bit の積を計算し、その後剰余を取る
        // これにより ulong の範囲を超える中間結果によるオーバーフローを防ぐ
        UInt128 product = (UInt128)a * b;
        return (ulong)(product % m);
    }

    /// <summary>
    /// (baseVal ^ exponent) mod modulus を計算します (繰り返し二乗法)。
    /// </summary>
    private static ulong ModPow(ulong baseVal, ulong exponent, ulong modulus)
    {
        if (modulus == 1) return 0; // mod 1 の結果は常に 0

        ulong result = 1;
        baseVal %= modulus; // 事前に基数を剰余で減らす

        while (exponent > 0)
        {
            // exponent の最下位ビットが 1 ならば、現在の baseVal を結果に乗算する
            if ((exponent & 1) == 1)
            {
                result = ModMultiply(result, baseVal, modulus);
            }

            // baseVal を二乗し、剰余を取る
            baseVal = ModMultiply(baseVal, baseVal, modulus);

            // exponent を右に 1 ビットシフト (exponent /= 2)
            exponent >>= 1;
        }

        return result;
    }
}

namespace CrimsonTeaLib
{
    public static class InputUtility
    {
        private static string[]? s_inputs;
        private static string? s_raw;
        private static int s_index = 0;

        private static void Init() => s_index = 0;
        public static int NextInt32 => int.Parse(s_inputs![s_index++]!);
        public static uint NextUInt32 => uint.Parse(s_inputs![s_index++]!);
        public static long NextInt64 => long.Parse(s_inputs![s_index++]!);
        public static BigInteger NextBigInt => BigInteger.Parse(s_inputs![s_index++]!);
        public static string NextString => s_inputs![s_index++];
        public static char NextChar => s_inputs![s_index++][0];
        public static int[] GetInt32Array() => s_inputs!.Select(int.Parse).ToArray();
        public static long[] GetInt64Array() => s_inputs!.Select(long.Parse).ToArray();
        public static string GetRawString() => s_raw!;
#if DEBUG
        private static TextReader? s_textReader;
        public static void SetSource(string path) => s_textReader = new StringReader(File.ReadAllText(path));
#endif
        public static bool InputNewLine()
        {
#if DEBUG
            if (s_textReader is TextReader sr)
            {
                Init();
                s_raw = sr.ReadLine()!;
                s_inputs = s_raw.Split(' ', StringSplitOptions.RemoveEmptyEntries);
                return true;
            }
#endif
            Init();
            s_raw = Console.ReadLine()!;
            s_inputs = s_raw.Split(' ', StringSplitOptions.RemoveEmptyEntries);
            return true;
        }
    }

    public static class CombinationFunction
    {
        public static void PartedRotate<T>(T[] a, int first1, int last1, int first2, int last2)
        {
            if (first1 == last1 || first2 == last2) return;
            int next = first2;
            while (first1 != next)
            {
                Swap(a, first1++, next++);
                if (first1 == last1) first1 = first2;
                if (next == last2)
                {
                    next = first2;
                }
                else if (first1 == first2)
                {
                    first2 = next;
                }
            }
        }

        public static bool NextCombinationImp<T>(T[] a, int first1, int last1, int first2, int last2) where T : IComparable<T>
        {
            if (first1 == last1 || first2 == last2) return false;
            int target = last1 - 1;
            int lastElem = last2 - 1;

            while (target != first1 && !(a[target].CompareTo(a[lastElem]) < 0)) target--;
            if (target == first1 && !(a[target].CompareTo(a[lastElem]) < 0))
            {
                PartedRotate(a, first1, last1, first2, last2);
                return false;
            }

            int next = a.AsSpan()[first2..].UpperBound(a[target]) + first2;
            Swap(a, target++, next++);
            PartedRotate(a, target, last1, next, last2);
            return true;
        }

        public static bool NextCombination<T>(T[] a, int first, int mid, int last) where T : IComparable<T>
            => NextCombinationImp(a, first, mid, mid, last);

        public static bool PrevCombination<T>(T[] a, int first, int mid, int last) where T : IComparable<T>
            => NextCombinationImp(a, mid, last, first, mid);

        public static void Swap<T>(T[] a, int i, int j) => (a[i], a[j]) = (a[j], a[i]);
    }

    public static class PermutationFunction
    {
        public static bool NextPermutation<T>(T[] a) where T : IComparable<T>
        {
            int n = a.Length;

            int i = n - 2;
            while (i >= 0 && a[i].CompareTo(a[i + 1]) >= 0) { i--; }

            if (i < 0) { return false; }

            int j = n - 1;
            while (a[j].CompareTo(a[i]) <= 0) { j--; }
            (a[i], a[j]) = (a[j], a[i]);

            Array.Reverse(a, i + 1, n - i - 1);

            return true;
        }
    }

    public readonly struct Output : IDisposable
    {
        private readonly StreamWriter _sw;
#if DEBUG
        public Output(string path)
        {
            var fs = new FileStream(path, FileMode.Create, FileAccess.Write);
            _sw = new StreamWriter(fs);
            Console.SetOut(_sw);
        }
#endif
        public Output(bool autoFlush)
        {
            _sw = new StreamWriter(Console.OpenStandardOutput()) { AutoFlush = autoFlush };
            Console.SetOut(_sw);
        }

        public void Dispose()
        {
            _sw.Dispose();
        }

        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        internal void Flush()
        {
            _sw.Flush();
        }
    }

    public static class ArrayExtensions
    {
        public static void Swap<T>(this T[] array, int i, int j) => (array[i], array[j]) = (array[j], array[i]);

        public static int LowerBound<T>(this T[] a, T target) where T : IComparable<T>
        {
            int ok = a.Length;
            int ng = -1;

            while (Math.Abs(ok - ng) > 1 && (ok + ng) / 2 is int mid)
                (ok, ng) = a[mid].CompareTo(target) >= 0 ? (mid, ng) : (ok, mid);

            return ok;
        }

        public static int UpperBound<T>(this T[] a, T target) where T : IComparable<T>
        {
            int ok = a.Length;
            int ng = -1;

            while (Math.Abs(ok - ng) > 1 && (ok + ng) / 2 is int mid)
                (ok, ng) = a[mid].CompareTo(target) > 0 ? (mid, ng) : (ok, mid);

            return ok;
        }

        public static int LowerBound<T>(this Span<T> a, T target) where T : IComparable<T>
        {
            int ok = a.Length;
            int ng = -1;

            while (Math.Abs(ok - ng) > 1 && (ok + ng) / 2 is int mid)
                (ok, ng) = a[mid].CompareTo(target) >= 0 ? (mid, ng) : (ok, mid);

            return ok;
        }

        public static int UpperBound<T>(this Span<T> a, T target) where T : IComparable<T>
        {
            int ok = a.Length;
            int ng = -1;

            while (Math.Abs(ok - ng) > 1 && (ok + ng) / 2 is int mid)
                (ok, ng) = a[mid].CompareTo(target) > 0 ? (mid, ng) : (ok, mid);

            return ok;
        }

        public struct IndexedEnumerable<T> : IEnumerable<(T item, int index)>
        {
            private readonly T[] _a;
            private readonly int _startIndex;

            public IndexedEnumerable(T[] a, int startIndex = 0)
            {
                _a = a;
                _startIndex = startIndex;
            }

            public readonly IndexedEnumerator<T> GetEnumerator() => new IndexedEnumerator<T>(_a, _startIndex);
            IEnumerator IEnumerable.GetEnumerator() => GetEnumerator();
            IEnumerator<(T item, int index)> IEnumerable<(T item, int index)>.GetEnumerator() => GetEnumerator();
        }

        public struct IndexedEnumerator<T> : IEnumerator<(T item, int index)>
        {
            public readonly (T item, int index) Current => (_a[_index], _index + _startIndex);

            private int _index;
            private int _startIndex;
            private T[] _a;

            public IndexedEnumerator(T[] a, int startIndex)
            {
                _index = -1;
                _a = a;
                _startIndex = startIndex;
            }

            public bool MoveNext() => ++_index < _a.Length;
            readonly object IEnumerator.Current => Current;
            public readonly void Dispose() { }
            public void Reset() => _index = -1;
        }

        /// <returns>(T value, int index)</returns>
        public static IndexedEnumerable<T> Enumerate<T>(this T[] arr, int startIndex = 0) => new IndexedEnumerable<T>(arr, startIndex);
    }

    public static class IEnumerableExtensions
    {
        public static IEnumerable<TSource> Log<TSource>(this IEnumerable<TSource> source)
        {
            Console.WriteLine(string.Join(' ', source));
            return source;
        }

        public static ScanEnumerable<TSource, TAccumulate> Scan<TSource, TAccumulate>(
            this IEnumerable<TSource> source,
            TAccumulate seed,
        Func<TAccumulate, TSource, TAccumulate> accumulator) where TSource : struct where TAccumulate : struct
        {
            return new ScanEnumerable<TSource, TAccumulate>(source, accumulator, seed);
        }

        public static IEnumerable<TAccumulate> ScanExSeed<TSource, TAccumulate>(
            this IEnumerable<TSource> source,
            TAccumulate seed,
            Func<TAccumulate, TSource, TAccumulate> accumulator)
        {
            var accumulation = new List<TAccumulate>();

            var current = seed;
            foreach (var item in source)
            {
                current = accumulator(current, item);
                accumulation.Add(current);
            }

            return accumulation;
        }

        public readonly struct ScanEnumerable<TSource, TAccumulate> : IEnumerable<TAccumulate> where TSource : struct where TAccumulate : struct
        {
            private readonly IEnumerable<TSource> _source;
            private readonly Func<TAccumulate, TSource, TAccumulate> _accumulator;
            private readonly TAccumulate _seed;

            public ScanEnumerable(IEnumerable<TSource> source, Func<TAccumulate, TSource, TAccumulate> accumulator, TAccumulate seed)
            {
                _source = source;
                _accumulator = accumulator;
                _seed = seed;
            }

            public readonly ScanEnumerator<TSource, TAccumulate> GetEnumerator() => new(_source, _accumulator, _seed);
            readonly IEnumerator IEnumerable.GetEnumerator() => GetEnumerator();
            readonly IEnumerator<TAccumulate> IEnumerable<TAccumulate>.GetEnumerator() => GetEnumerator();
        }

        public struct ScanEnumerator<TSource, TAccumulate> : IEnumerator<TAccumulate> where TSource : struct where TAccumulate : struct
        {
            private readonly Func<TAccumulate, TSource, TAccumulate> _accumulator;

            private readonly IEnumerator<TSource> _enumerator;
            private TAccumulate _current;

            private bool _secondOrLaterElement = false;

            public ScanEnumerator(IEnumerable<TSource> source, Func<TAccumulate, TSource, TAccumulate> accumulator, TAccumulate seed)
            {
                _enumerator = source.GetEnumerator();
                _accumulator = accumulator;

                _current = seed;
            }

            public readonly TAccumulate Current => _current;
            readonly object IEnumerator.Current => Current;

            public readonly void Dispose() { }

            public bool MoveNext()
            {
                if (_secondOrLaterElement)
                {
                    if (_enumerator.MoveNext())
                    {
                        _current = _accumulator(_current, _enumerator.Current);
                        return true;
                    }
                    return false;
                }
                else
                {
                    _secondOrLaterElement = true;
                    return true;
                }
            }

            public void Reset()
            {
                throw new NotSupportedException();
            }
        }

        public static IEnumerable<TSource> Scan<TSource>(
            this IEnumerable<TSource> source,
            Func<TSource, TSource, TSource> accumulator)
        {
            if (source is null)
                throw new ArgumentNullException(paramName: nameof(source));

            if (accumulator is null)
                throw new ArgumentNullException(paramName: nameof(accumulator));

            var accumulation = new List<TSource>();

            if (source.Any() is false)
            {
                return accumulation;
            }

            var current = source.First();
            accumulation.Add(current);

            foreach (var item in source.Skip(1))
            {
                current = accumulator(current, item);
                accumulation.Add(current);
            }

            return accumulation;
        }

        public static CombinationEnumerable<T> Combination<T>(this T[] a, int k) where T : IComparable<T>
            => new(a, k);

        public readonly struct CombinationEnumerable<T> where T : IComparable<T>
        {
            private readonly T[] _a;
            private readonly int _k;

            public CombinationEnumerable(T[] a, int k)
            {
                _a = a;
                _k = k;
            }

            public readonly CombinationEnumerator<T> GetEnumerator() => new(_a, _k);
        }

        public struct CombinationEnumerator<T> : IEnumerator<ReadOnlyMemory<T>> where T : IComparable<T>
        {
            private readonly int _k;
            private readonly T[] _a;
            private readonly int _n;
            private bool _secondOrLaterElement = false;

            public CombinationEnumerator(T[] a, int k)
            {
                _a = a;
                _n = a.Length;
                _k = k;
            }

            public readonly ReadOnlyMemory<T> Current => _a.AsMemory()[.._k];
            readonly object IEnumerator.Current => Current;

            public readonly void Dispose() { }

            public bool MoveNext()
            {
                if (_secondOrLaterElement)
                {
                    return CombinationFunction.NextCombination(_a, 0, _k, _n);
                }
                else
                {
                    _secondOrLaterElement = true;
                    return true;
                }
            }

            public void Reset()
            {
                throw new NotSupportedException();
            }
        }

        public static PermutationEnumerable<T> Permutation<T>(this T[] a) where T : IComparable<T> => new(a);

        public readonly struct PermutationEnumerable<T> where T : IComparable<T>
        {
            private readonly T[] _a;
            public PermutationEnumerable(T[] a) => _a = a;
            public readonly PermutationEnumerator<T> GetEnumerator() => new(_a);
        }

        public struct PermutationEnumerator<T> : IEnumerator<ReadOnlyMemory<T>> where T : IComparable<T>
        {
            private readonly T[] _a;
            private readonly int _n;
            private bool _secondOrLaterElement = false;

            public PermutationEnumerator(T[] a)
            {
                _a = a;
                _n = a.Length;
            }
            public readonly ReadOnlyMemory<T> Current => _a.AsMemory()[..];
            readonly object IEnumerator.Current => Current;

            public readonly void Dispose() { }

            public bool MoveNext()
            {
                if (_secondOrLaterElement)
                {
                    return PermutationFunction.NextPermutation(_a);
                }
                else
                {
                    _secondOrLaterElement = true;
                    return true;
                }
            }

            public void Reset() => throw new NotSupportedException();
        }
    }

    public class MyMath
    {
        public static long Pow(long x, int y) => Enumerable.Repeat(x, y).Aggregate(1L, (acc, x) => acc * x);
        public static long Pow10(int y) => Pow(10, y);
        public static long Pow2(int y) => Pow(2, y);
    }

    public class MyMathBigInteger
    {
        public static BigInteger Pow(long x, int y) => Enumerable.Repeat(x, y).Aggregate(new BigInteger(1), (acc, x) => acc * x);
        public static BigInteger Pow10(int y) => Pow(10, y);
        public static BigInteger Pow2(int y) => Pow(2, y);
    }
}

yukicoder

結果

コンパイルメッセージ

ソースコード