#https://scrapbox.io/wandsbox/%E7%B4%A0%E6%95%B0%E3%81%AE%E4%B8%80%E8%88%AC%E9%A0%85%E3%82%92%E6%B1%82%E3%82%81%E3%82%88

import math
import random
from sympy import prime, primepi, isprime


#Naive: 遅いがシンプル

#エラトステネスの篩（Sieve of Eratosthenes）: 小さい範囲（例: 10⁶以下）の素数を列挙したい
#最速かつ実装が簡単。
#時間計算量: 𝑂(𝑛log log 𝑛)
#実装の工夫（ビットベース、奇数のみ処理）で非常に高速になる。

#AKS 素数判定法（理論的に決定的）: 非常に大きな素数（数百桁〜数千桁）の判定がしたい
#決定的素数判定法だが、現実的にはまだ遅い。

#Miller–Rabin: 極めて大きい数に対して有効（確率的な保証）
#指数的に高速かつ十分な確率精度を持つ。

#Generator: 必要なだけ素数を取り出せる（無限列）


def is_prime_naive(n):
    if n < 2:
        return False
    for i in range(2, int(n**0.5)+1):
        if n % i == 0:
            return False
    return True

def primes_naive(limit):
    return [n for n in range(2, limit) if is_prime_naive(n)]

def sieve_eratosthenes(limit):
    sieve = [True] * limit
    sieve[0:2] = [False, False]  # 0と1は素数ではない
    for i in range(2, int(limit**0.5)+1):
        if sieve[i]:
            sieve[i*i:limit:i] = [False] * len(range(i*i, limit, i))
    return [i for i, is_prime in enumerate(sieve) if is_prime]

def is_prime_miller_rabin(n, k=5):
    if n <= 1:
        return False
    if n <= 3:
        return True
    if n % 2 == 0:
        return False

    # n - 1 = 2^r * d の形に変形
    r, d = 0, n - 1
    while d % 2 == 0:
        d //= 2
        r += 1

    for _ in range(k):  # kは繰り返し数（信頼度）
        a = random.randrange(2, n - 2)
        x = pow(a, d, n)

        if x in (1, n - 1):
            continue

        for _ in range(r - 1):
            x = pow(x, 2, n)
            if x == n - 1:
                break
        else:
            return False
    return True

def primes_miller_rabin(limit):
    return [n for n in range(2, limit) if is_prime_miller_rabin(n)]

def gen_primes():
    primes = []
    n = 2
    while True:
        for p in primes:
            if n % p == 0:
                break
        else:
            primes.append(n)
            yield n
        n += 1

def first_n_primes(n):
    result = []
    for i, p in enumerate(gen_primes()):
        if i >= n:
            break
        result.append(p)
    return result

def euler_prime_polynomial(n_max):
    primes = []
    for n in range(n_max):
        val = n * n + n + 41
        if is_prime_naive(val):  # 上のnaive判定を使う
            primes.append(val)
    return primes

def is_prime(num):
    """
    指定された数が素数かどうかを判定する関数
    (試行除算を使用)
    """
    if num < 2:
        return False
    for i in range(2, int(math.sqrt(num)) + 1):
        if num % i == 0:
            return False
    return True

def count_primes_up_to(limit):
    """
    指定された上限までの素数の個数を数える関数
    """
    count = 0
    for i in range(2, limit + 1):
        if is_prime(i):
            count += 1
    return count

def get_nth_prime_by_sum_method(n, max_search_limit):
    """
    Args:
        n (int): 求めたい素数の順番 (例: 10で10番目の素数)
        max_search_limit (int): 探索する最大の自然数 (f(n)に相当)

    Returns:
        int: n番目の素数。見つからなかった場合はNone。
    """
    if n <= 0:
        raise ValueError("nは正の整数である必要があります。")

    total_sum = 0
    for m in range(1, max_search_limit + 1):
        is_m_prime = is_prime(m)
        num_primes_below_m = count_primes_up_to(m)

        # mが素数 かつ m以下の素数がn個であるならばmを加算
        if is_m_prime and num_primes_below_m == n:
            total_sum += m
            # n番目の素数が見つかった時点でループを終了
            # このアルゴリズムの性質上、n番目の素数 'のみ' がこの条件を満たすため、
            # 最初に見つかった時点でそれがn番目の素数であり、それ以降のmは条件を満たさない。
            return total_sum
    return None # 指定された範囲内でn番目の素数が見つからなかった場合

# n番目の素数がおよそ n*ln(n) くらいになることを考慮して、max_search_limitを設定
def estimate_upper_bound(n):
    if n < 6: # 小さいnには固定値を設定
        return 15
    return int(n * math.log(n) * 1.5) + 10 # 経験的な係数


if __name__ == "__main__":
    print("Naive:", primes_naive(50))
    print("Eratosthenes:", sieve_eratosthenes(50))
    print("Miller-Rabin:", primes_miller_rabin(50))
    print("Generator:", first_n_primes(20))
    print("euler_prime_polynomial:", euler_prime_polynomial(20))

    """
# n 以下の素数をすべてリストで取得
    primes = primesieve.primes(100)
    print("Primes up to 100:", primes)

# n から m の範囲の素数を取得
    primes_range = primesieve.primes(100, 200)
    print("Primes from 100 to 200:", primes_range)

# n 以下の素数の個数（π(n)）
    count = primesieve.count_primes(0, 1_000_000)
    print("π(1,000,000) =", count)

# n 番目の素数（1番目は2）
    nth = primesieve.nth_prime(10000)
    print("10000th prime =", nth)
    """


# n 番目の素数（1番目は2）
    p = prime(10000)
    print("10000th prime =", p)

# n 以下の素数の個数（π(n)）
    count = primepi(1_000_000)
    print("π(1,000,000) =", count)

# 素数判定
    print("Is 104729 prime?", isprime(104729))

    n_value = 10
    limit = estimate_upper_bound(n_value)
    result = get_nth_prime_by_sum_method(n_value, limit)
    print(f"{n_value}番目の素数は {result} です。")

    n_value = 100
    limit = estimate_upper_bound(n_value)
    result = get_nth_prime_by_sum_method(n_value, limit)
    print(f"{n_value}番目の素数は {result} です。")

    n_value = 200
    limit = estimate_upper_bound(n_value)
    result = get_nth_prime_by_sum_method(n_value, limit)
    print(f"{n_value}番目の素数は {result} です。")