// Copyright (C) 2026 Kiyotsugu Arai
// SPDX-License-Identifier: LGPL-3.0-or-later

// IntSqrt.hpp
// 整数平方根・n乗根と完全冪判定
//
// アルゴリズム:
// - sqrt: mpn レベル Newton 法 (小さい値は double 精度)
// - sqrtRem: 平方根と余りを同時に計算
// - nthRoot: Precision Doubling Newton 法によるn乗根計算
// - nthRootRem: n乗根と余りを同時に計算
// - isSquare: GMP 式 3 層フィルタ (mod 256 + mod_34lsub1 + sqrtRem)
//
// 参考文献:
// - Brent & Zimmermann, "Modern Computer Arithmetic", §1.5.2
// - GMP: mpz_sqrt, mpz_sqrtrem, mpz_root, mpz_rootrem, mpz_perfect_square_p

#ifndef CALX_INT_SQRT_HPP
#define CALX_INT_SQRT_HPP

#include "IntBase.hpp"
#include <utility>

namespace calx {

class IntSqrt {
public:
    // 整数平方根 (floor(√value))
    // 特殊状態: NaN → NaN, 負 → NaN (NegativeSqrt), +∞ → +∞, 0 → 0
    // アルゴリズム: サイズに応じて自動選択
    //   - MSB < 52: double 精度で計算 (正確)
    //   - それ以上: mpn レベル Newton 法
    static Int sqrt(const Int& value);

    // 平方根と余りを返す: value = sqrt² + remainder
    // remainder は常に 0 ≤ remainder < 2*sqrt + 1 を満たす
    // GMP 互換: mpz_sqrtrem
    static Int sqrtRem(const Int& value, Int& remainder);

    // 完全平方数判定
    // pSqrt != nullptr の場合、平方根を格納する
    // 高速化: GMP 式 3 層フィルタ (mod 256, mod_34lsub1, sqrtRem)
    // GMP 互換: mpz_perfect_square_p
    static bool isSquare(const Int& value, Int* pSqrt = nullptr);

    // n乗根 (floor(value^(1/n)))
    // 特殊状態: NaN → NaN, n≤0 → NaN, 負でnが偶数 → NaN, +∞ → +∞, 0 → 0
    // アルゴリズム: Precision Doubling Newton + remainder 追跡
    //   x_{k+1} = ((n-1)*x_k + value/x_k^(n-1)) / n
    // GMP 互換: mpz_root
    static Int nthRoot(const Int& value, uint32_t n);

    // n乗根と余りを返す: value = root^n + remainder
    // remainder は常に 0 ≤ remainder を満たす
    // GMP 互換: mpz_rootrem
    static Int nthRootRem(const Int& value, uint32_t n, Int& remainder);

    // 完全冪判定: value = b^k (k >= 2) となる k が存在するか
    // 存在すれば true を返し、pBase != nullptr なら b を格納
    // アルゴリズム: k = 2..bitLength を試行、nthRoot で検証
    static bool isPerfectPower(const Int& value, Int* pBase = nullptr, uint32_t* pExp = nullptr);

private:
    // double 精度で平方根計算 (MSB < 52 ビット用)
    static Int sqrt_small(const Int& value);
};

} // namespace calx

#endif // CALX_INT_SQRT_HPP