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

// Float.hpp
// 多倍長浮動小数点数クラス
//
// このファイルでは、多倍長整数（Int）を仮数部として使用する
// 高精度浮動小数点演算のためのクラスを提供します。
//
// 主な機能:
// - 任意精度の浮動小数点数
// - 標準的な算術演算と数学関数
// - IEEE 754に準拠した丸めモード
// - 高度な数値演算機能

#ifndef CALX_FLOAT_HPP
#define CALX_FLOAT_HPP

#include <math/core/mp/Int.hpp>
#include <math/concepts/algebraic_concepts.hpp>
#include <math/core/common.hpp>
#include <math/core/traits.hpp>

#include <string>
#include <string_view>
#include <iosfwd>
#include <compare>  // C++20 三方比較演算子
#include <span>
#include <climits>      // INT_MAX
#include <type_traits>  // std::is_integral_v
#include <utility>      // std::pair

namespace calx {

    // 精度伝播方式の選択
    // MIN_PROPAGATION: lib++20 方式 (min) — 段階的精度、計算量節約
    // MAX_PROPAGATION: lib++23 方式 (max/merge) — 精度維持、安全性重視
    enum class PrecisionPolicy { MIN_PROPAGATION, MAX_PROPAGATION };

#ifdef CALX_PRECISION_POLICY_MIN
    inline constexpr PrecisionPolicy PRECISION_POLICY = PrecisionPolicy::MIN_PROPAGATION;
#else
    inline constexpr PrecisionPolicy PRECISION_POLICY = PrecisionPolicy::MAX_PROPAGATION;
#endif

    // Float クラスの前方宣言
    class Float;

    // 数学関数の前方宣言
    [[nodiscard]] Float exp(const Float& x, int precision);
    [[nodiscard]] Float exp(Float&& x, int precision);
    [[nodiscard]] Float log(const Float& x, int precision);
    [[nodiscard]] Float log(Float&& x, int precision);
    [[nodiscard]] Float logUi(unsigned long long n, int precision);
    [[nodiscard]] Float sin(const Float& x, int precision);
    [[nodiscard]] Float sin(Float&& x, int precision);
    [[nodiscard]] Float cos(const Float& x, int precision);
    [[nodiscard]] Float cos(Float&& x, int precision);
    [[nodiscard]] Float tan(const Float& x, int precision);
    [[nodiscard]] Float tan(Float&& x, int precision);
    [[nodiscard]] Float sqr(const Float& x, int precision);
    [[nodiscard]] Float sqr(Float&& x, int precision);
    [[nodiscard]] Float sqrt(const Float& x, int precision);
    [[nodiscard]] Float sqrt(Float&& x, int precision);
    [[nodiscard]] Float sinh(const Float& x, int precision);
    [[nodiscard]] Float sinh(Float&& x, int precision);
    [[nodiscard]] Float cosh(const Float& x, int precision);
    [[nodiscard]] Float cosh(Float&& x, int precision);
    [[nodiscard]] Float tanh(const Float& x, int precision);
    [[nodiscard]] Float tanh(Float&& x, int precision);
    [[nodiscard]] Float pow(const Float& x, const Float& y, int precision);
    [[nodiscard]] Float pow(Float&& x, Float&& y, int precision);
    [[nodiscard]] Float pow(const Float& x, int n, int precision);
    [[nodiscard]] Float pow(Float&& x, int n, int precision);
    [[nodiscard]] Float abs(const Float& value);
    [[nodiscard]] Float abs(Float&& value);
    [[nodiscard]] Float ldexp(const Float& value, int exp);
    [[nodiscard]] Float ldexp(Float&& value, int exp);
    [[nodiscard]] Float frexp(const Float& value, int* exp);
    [[nodiscard]] Float modf(const Float& value, Float* iptr);

    // 丸めモード（IEEE 754に準拠）
    enum class RoundingMode {
        ToNearest,      // 最も近い値に丸め（デフォルト）
        TowardZero,     // ゼロ方向に丸め（切り捨て）
        TowardPositive, // 正の無限大方向に丸め（切り上げ）
        TowardNegative, // 負の無限大方向に丸め（切り下げ）
        AwayFromZero    // ゼロから遠ざかる方向に丸め
    };

    // IEEE 754 例外フラグ
    enum FloatException : unsigned {
        FE_NONE      = 0,
        FE_INEXACT   = 1u << 0,   // 不正確 (丸めが発生)
        FE_UNDERFLOW = 1u << 1,   // アンダーフロー
        FE_OVERFLOW  = 1u << 2,   // オーバーフロー
        FE_DIVBYZERO = 1u << 3,   // ゼロ除算
        FE_INVALID   = 1u << 4,   // 無効演算 (NaN 生成)
        FE_ALL       = 0x1Fu
    };

    /**
     * @brief 多倍長浮動小数点数クラス
     *
     * この実装では、多倍長整数（Int）を仮数部として使用し、
     * 別途指数部を保持することで高精度な浮動小数点演算を可能にします。
     * IEEE 754の精度と丸めモードの概念に準拠しています。
     */
    class Float {
        // 友達クラスの宣言
        friend struct numeric_traits<Float>;
        friend class FloatOps;

        // 演算子のフレンド宣言
        friend Float operator+(const Float& value);  // 単項プラス演算子
        friend Float operator-(const Float& value);  // 単項マイナス演算子
        friend Float operator+(const Float& lhs, const Float& rhs);
        friend Float operator+(Float&& lhs, const Float& rhs);
        friend Float operator+(const Float& lhs, Float&& rhs);
        friend Float operator+(Float&& lhs, Float&& rhs);
        friend Float operator-(const Float& lhs, const Float& rhs);
        friend Float operator-(Float&& lhs, const Float& rhs);
        friend Float operator-(const Float& lhs, Float&& rhs);
        friend Float operator-(Float&& lhs, Float&& rhs);
        friend Float operator*(const Float& lhs, const Float& rhs);
        friend Float operator*(Float&& lhs, const Float& rhs);
        friend Float operator*(const Float& lhs, Float&& rhs);
        friend Float operator*(Float&& lhs, Float&& rhs);
        friend Float operator/(const Float& lhs, const Float& rhs);
        friend Float operator/(Float&& lhs, const Float& rhs);
        friend Float mulScalarF(const Float& lhs, uint64_t rhs);
        friend Float mulScalarF(Float&& lhs, uint64_t rhs);
        friend Float mulScalarF(const Float& lhs, int64_t rhs);
        friend Float mulScalarF(Float&& lhs, int64_t rhs);
        friend Float divScalarF(const Float& lhs, uint64_t rhs);
        friend Float divScalarF(Float&& lhs, uint64_t rhs);
        friend Float divScalarF(const Float& lhs, int64_t rhs);
        friend Float divScalarF(Float&& lhs, int64_t rhs);
        friend std::partial_ordering operator<=>(const Float& lhs, const Float& rhs);
        friend bool operator==(const Float& lhs, const Float& rhs);

        // 数学関数のフレンド宣言
        friend Float abs(const Float& value);
        friend Float abs(Float&& value);
        friend Float exp(const Float& x, int precision);
        friend Float log(const Float& x, int precision);
        friend Float sin(const Float& x, int precision);
        friend Float cos(const Float& x, int precision);
        friend Float tan(const Float& x, int precision);
        friend Float sqr(const Float& x, int precision);
        friend Float sqrt(const Float& x, int precision);
        friend Float importBinaryFloat(std::span<const uint8_t> data);
        friend Float sinh(const Float& x, int precision);
        friend Float cosh(const Float& x, int precision);
        friend Float tanh(const Float& x, int precision);
        friend Float pow(const Float& x, const Float& y, int precision);
        friend Float pow(const Float& x, int n, int precision);
        friend Float ldexp(const Float& value, int exp);
        friend Float ldexp(Float&& value, int exp);
        friend Float frexp(const Float& value, int* exp);

        // ユーティリティ関数のフレンド宣言
        friend void swap(Float& a, Float& b) noexcept;
        friend Float copySign(const Float& x, const Float& y);
        friend Float copySign(Float&& x, const Float& y);

    private:
        Int mantissa_;       // 仮数部（常に非負。符号は is_negative_ が管理）
        int64_t exponent_;   // 指数部
        bool is_negative_;   // 符号フラグ（Int の符号とは独立）
            // mantissa_ は常に m_sign >= 0 に保つ。符号を別に持つ理由:
            // (1) mpn 関数が符号なし配列を前提とし、符号付き分岐が不要になる
            // (2) IEEE 754 の -0.0 を表現するため (Int では zero の符号は常に 0)

        // 特殊値フラグ
        bool is_infinity_;   // 無限大フラグ
        bool is_nan_;        // NaNフラグ

        // 精度管理フィールド
        // effective_bits_: 仮数に実際に含まれている信頼できるビット数
        //   - 乗除算: min(lhs, rhs)
        //   - 加算: min(lhs, rhs) (桁落ちなし)
        //   - 減算: min(lhs, rhs) - lost_bits (桁落ち考慮)
        //   - 整数由来 (exact): INT_MAX
        //   - double 由来: 53
        // requested_bits_: この値に求められる目標精度
        //   - 演算結果: max(lhs, rhs)
        //   - setPrecision() で明示的に設定可能
        //   - 整数由来 (exact): INT_MAX
        int effective_bits_; // 有効ビット数 (信頼できる精度)
        int requested_bits_; // 要求ビット数 (目標精度)

        // 精度とモード設定
        // inline: ヘッダ定義で ODR 準拠。DLL 境界を越えても呼び出し側の TU に変数が存在する。
        inline static thread_local int default_precision_ = 53;  // デフォルトはdoubleと同等の精度
        inline static thread_local RoundingMode rounding_mode_ = RoundingMode::ToNearest;
        static thread_local int64_t emin_;                 // 最小指数 (IEEE 754 emin)
        static thread_local int64_t emax_;                 // 最大指数 (IEEE 754 emax)
        static thread_local unsigned exception_flags_;     // IEEE 754 例外フラグ

        // 指数部のオーバーフロー防御定数
        // この範囲を超えると ±∞ または 0 に丸められる
        static constexpr int64_t EXPONENT_MAX = (1LL << 60);
        static constexpr int64_t EXPONENT_MIN = -(1LL << 60);

        // 正規化関数（仮数部を調整し指数部を更新）
        void normalize();

        // 指数部の範囲チェック（オーバーフロー → ±∞、アンダーフロー → 0）
        void checkExponentBounds();

    public:
        // 型定義（標準コンテナ互換）
        using value_type = Float;
        using size_type = std::size_t;
        using difference_type = std::ptrdiff_t;
        using reference = Float&;
        using const_reference = const Float&;

        // コンストラクタとデストラクタ

        /**
         * @brief デフォルトコンストラクタ（ゼロに初期化）
         */
        Float();

        /**
         * @brief 整数からのコンストラクタ
         * @param value 初期値
         */
        explicit Float(int value);

        /**
         * @brief 64ビット整数からのコンストラクタ
         * @param value 初期値
         */
        explicit Float(int64_t value);

        /**
         * @brief 整数仮数部と指数部を指定するコンストラクタ
         * @param mantissa 整数仮数部
         * @param exponent 指数部
         * @param is_negative 符号フラグ
         */
        explicit Float(int64_t mantissa, int64_t exponent, bool is_negative = false);

        /**
         * @brief 多倍長整数からのコンストラクタ
         * @param value 初期値
         */
        explicit Float(const Int& value);

        /**
         * @brief 浮動小数点数からのコンストラクタ
         * @param value 初期値
         */
        explicit Float(double value);

        /**
         * @brief 文字列からのコンストラクタ
         * @param str 数値を表す文字列
         */
        explicit Float(std::string_view str);

        /**
         * @brief 仮数部と指数部を指定するコンストラクタ
         * @param mantissa 仮数部
         * @param exponent 指数部
         * @param is_negative 符号フラグ
         */
        Float(const Int& mantissa, int64_t exponent, bool is_negative = false);
        Float(Int&& mantissa, int64_t exponent, bool is_negative = false);

        /**
         * @brief コピーコンストラクタ
         */
        Float(const Float& other) = default;

        /**
         * @brief ムーブコンストラクタ
         */
        Float(Float&& other) noexcept = default;

        /**
         * @brief コピー代入演算子
         */
        Float& operator=(const Float& other) = default;

        /**
         * @brief ムーブ代入演算子
         */
        Float& operator=(Float&& other) noexcept = default;

        /**
         * @brief 整数値からの代入演算子
         */
        Float& operator=(int64_t value);

        /**
         * @brief 浮動小数点数からの代入演算子
         */
        Float& operator=(double value);

        /**
         * @brief 加算代入演算子
         */
        Float& operator+=(const Float& rhs);

        /**
         * @brief 減算代入演算子
         */
        Float& operator-=(const Float& rhs);

        /**
         * @brief 乗算代入演算子
         */
        Float& operator*=(const Float& rhs);

        /**
         * @brief 除算代入演算子
         */
        Float& operator/=(const Float& rhs);

        /**
         * @brief デストラクタ
         */
        ~Float() = default;

        // 基本アクセス関数

        /**
         * @brief 仮数部を取得
         * @return 仮数部の参照
         */
        [[nodiscard]] const Int& mantissa() const { return mantissa_; }

        /**
         * @brief 指数部を取得
         * @return 指数部の値
         */
        [[nodiscard]] int64_t exponent() const { return exponent_; }

        /**
         * @brief 符号を取得
         * @return 負の数ならtrue
         */
        [[nodiscard]] bool isNegative() const { return is_negative_; }

        /**
         * @brief ゼロかどうかをチェック
         * @return ゼロならtrue
         */
        [[nodiscard]] bool isZero() const;

        /**
         * @brief 無限大かどうかをチェック
         * @return 無限大ならtrue
         */
        [[nodiscard]] bool isInfinity() const { return is_infinity_; }

        /**
         * @brief NaNかどうかをチェック
         * @return NaNならtrue
         */
        [[nodiscard]] bool isNaN() const { return is_nan_; }

        /**
         * @brief 正確な値（整数由来）かどうかをチェック
         * @return 正確な値ならtrue（演算時は相手の精度に適応する）
         */
        [[nodiscard]] bool isExact() const { return effective_bits_ >= INT_MAX && requested_bits_ >= INT_MAX; }

        /**
         * @brief 有効ビット数を取得
         * @return 有効ビット数
         */
        [[nodiscard]] int effectiveBits() const { return effective_bits_; }

        /// 有効ビット数を直接設定する (数学関数が入力精度を反映する際に使用)
        void setEffectiveBits(int bits) { effective_bits_ = bits; }

        /**
         * @brief 要求ビット数を取得
         * @return 要求ビット数
         */
        [[nodiscard]] int requestedBits() const { return requested_bits_; }

        /**
         * @brief 有効桁数を取得
         * @return 有効桁数
         */
        [[nodiscard]] int precision() const;

        /**
         * @brief 有効桁数を設定
         * @param precision 新しい有効桁数
         * @return 自身への参照
         */
        Float& setPrecision(int precision);

        /**
         * @brief 数学関数の結果に精度フィールドを設定
         * @param precision 計算精度 (10進桁数)
         * @return 自身への参照
         *
         * setPrecision(precision) を呼び出した後、effective_bits_ を
         * precisionToBits(precision) に設定する。数学関数が結果を返す直前に使用。
         */
        Float& setResultPrecision(int precision);

        /**
         * @brief ワード単位の高速近似切り詰め (中間計算用)
         * @param precision 目標精度 (10進桁数)
         *
         * setPrecision と異なり、丸め処理を行わずワード単位で下位を切り捨てる。
         * 最大63ビットの余剰を保持するが、ガードビットで吸収される。
         * ループ内の中間値の精度制御に使用し、最終結果には setPrecision を使う。
         */
        void truncateToApprox(int precision);

        // 精度とビット数の変換ユーティリティ
        [[nodiscard]] static int precisionToBits(int precision);
        [[nodiscard]] static int bitsToPrecision(int bits);

        /// 入力の requestedBits から適切な10進精度を決定するヘルパー
        /// 特殊値 (infinity, NaN, zero) や exact 値は defaultPrecision にフォールバック
        [[nodiscard]] static int requestedPrecision(const Float& x) {
            int rb = x.requestedBits();
            int dp = defaultPrecision();
            return (rb <= 0 || rb >= INT_MAX)
                   ? dp
                   : std::max(bitsToPrecision(rb), dp);
        }

        /**
         * @brief デフォルトの有効桁数を取得
         * @return デフォルト有効桁数
         */
        [[nodiscard]] static int defaultPrecision() { return default_precision_; }

        /**
         * @brief デフォルトの有効桁数を設定
         * @param precision 新しいデフォルト有効桁数
         */
        static void setDefaultPrecision(int precision) { default_precision_ = precision; }

        /**
         * @brief 丸めモードを取得
         * @return 現在の丸めモード
         */
        [[nodiscard]] static RoundingMode roundingMode() { return rounding_mode_; }

        /**
         * @brief 丸めモードを設定
         * @param mode 新しい丸めモード
         */
        static void setRoundingMode(RoundingMode mode) { rounding_mode_ = mode; }

        // ── 指数範囲制御 (IEEE 754 emin/emax) ──

        /// 最小指数を取得
        [[nodiscard]] static int64_t emin() { return emin_; }
        /// 最大指数を取得
        [[nodiscard]] static int64_t emax() { return emax_; }
        /// 最小指数を設定
        static void setEmin(int64_t e) { emin_ = e; }
        /// 最大指数を設定
        static void setEmax(int64_t e) { emax_ = e; }
        /// 指数範囲をデフォルト (±2^60) にリセット
        static void resetExponentRange() { emin_ = EXPONENT_MIN; emax_ = EXPONENT_MAX; }

        /// サブノーマルエミュレーション: 指数が emin を下回る場合に仮数部を右シフト
        void subnormalize(RoundingMode mode = rounding_mode_);

        // ── IEEE 754 例外フラグ管理 ──

        /// 例外フラグを取得
        [[nodiscard]] static unsigned getExceptionFlags() { return exception_flags_; }
        /// 指定した例外フラグをセット
        static void raiseException(unsigned flags) { exception_flags_ |= (flags & FE_ALL); }
        /// 指定した例外フラグをクリア
        static void clearExceptionFlags(unsigned flags = FE_ALL) { exception_flags_ &= ~(flags & FE_ALL); }
        /// 指定した例外が発生しているかテスト
        [[nodiscard]] static bool testException(unsigned flags) { return (exception_flags_ & flags) != 0; }

        /**
         * @brief 文字列に変換
         * @param precision 出力精度（デフォルトは現在の精度）
         * @return 文字列表現
         */
        [[nodiscard]] std::string toString(int precision = -1) const;

        /// N 進数表示 (例: 2進 "11.0101", 16進 "FF.A5")
        /// @param base      基数 (2〜36)
        /// @param fracDigits 小数部の桁数 (-1 で仮数部全ビットから自動決定)
        [[nodiscard]] std::string toString(int base, int fracDigits) const;

        /**
         * @brief 10進数の文字列に変換
         * @param precision 出力精度（デフォルトは現在の精度）
         * @return 10進数文字列
         */
        [[nodiscard]] std::string toDecimalString(int precision = -1) const;

        /**
         * @brief 科学表記の文字列に変換
         * @param precision 出力精度（デフォルトは現在の精度）
         * @return 科学表記文字列
         */
        [[nodiscard]] std::string toScientificString(int precision = -1) const;

        /**
         * @brief 浮動小数点数に変換（精度の損失の可能性あり）
         * @return double値
         */
        [[nodiscard]] double toDouble() const;

        // 特殊値の生成

        /**
         * @brief 正の無限大を生成
         * @return 正の無限大
         */
        [[nodiscard]] static Float positiveInfinity();

        /**
         * @brief 負の無限大を生成
         * @return 負の無限大
         */
        [[nodiscard]] static Float negativeInfinity();

        /**
         * @brief NaNを生成
         * @return NaN
         */
        [[nodiscard]] static Float nan();

        /**
         * @brief 機械イプシロンを生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 機械イプシロン
         */
        [[nodiscard]] static Float epsilon(int precision = default_precision_);

        /**
         * @brief 精度を指定してゼロを生成
         * @param precision 精度
         * @return 指定精度のゼロ
         */
        [[nodiscard]] static Float zero(int precision = default_precision_);

        /**
         * @brief 精度を指定して1を生成
         * @param precision 精度
         * @return 指定精度の1
         */
        [[nodiscard]] static Float one(int precision = default_precision_);

        /**
         * @brief 円周率πを生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度のπ
         */
        [[nodiscard]] static Float pi(int precision = default_precision_);

        /**
         * @brief 自然対数の底eを生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度のe
         */
        [[nodiscard]] static Float e(int precision = default_precision_);

        /**
         * @brief log(2)を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度のlog(2)
         */
        [[nodiscard]] static Float log2(int precision = default_precision_);

        /**
         * @brief log(10)を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度のlog(10)
         */
        [[nodiscard]] static Float log10(int precision = default_precision_);

        /**
         * @brief オイラー・マスケローニ定数 γ を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度の γ ≈ 0.5772156649...
         */
        [[nodiscard]] static Float euler(int precision = default_precision_);

        /**
         * @brief カタラン定数 G を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度の G ≈ 0.9159655941...
         */
        [[nodiscard]] static Float catalan(int precision = default_precision_);

        /**
         * @brief √2 を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度の √2 ≈ 1.4142135623...
         */
        [[nodiscard]] static Float sqrt2(int precision = default_precision_);

        /**
         * @brief レムニスケート定数 ω を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度の ω ≈ 2.6220575542...
         *
         * AGM を使用: ω = 4·(a+b) / (t·√2)
         * ここで a,b,t は AGM(1, 1/√2) の収束値
         */
        [[nodiscard]] static Float lemniscate(int precision = default_precision_);

        /**
         * @brief Γ(1/4) を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度の Γ(1/4) ≈ 3.6256099082...
         *
         * AGM を使用: Γ(1/4) = √(ω·√(2π))
         * = √(lemniscate · √(2π))
         */
        [[nodiscard]] static Float gamma14(int precision = default_precision_);

        /**
         * @brief Apéry 定数 ζ(3) を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度の ζ(3) ≈ 1.2020569031...
         *
         * Amdeberhan-Zeilberger 第2公式 (Binary Splitting)
         */
        [[nodiscard]] static Float zeta3(int precision = default_precision_);

        /**
         * @brief ζ(5) を生成
         * @param precision 精度（デフォルトは現在の精度）
         * @return 指定精度の ζ(5) ≈ 1.0369277551...
         *
         * Broadhurst BBP 公式 (Binary Splitting)
         */
        [[nodiscard]] static Float zeta5(int precision = default_precision_);

        // --- π 派生定数 ---
        [[nodiscard]] static Float half_pi(int precision = default_precision_);       // π/2
        [[nodiscard]] static Float quarter_pi(int precision = default_precision_);    // π/4
        [[nodiscard]] static Float two_pi(int precision = default_precision_);        // 2π (τ)
        [[nodiscard]] static Float inv_pi(int precision = default_precision_);        // 1/π
        [[nodiscard]] static Float two_inv_pi(int precision = default_precision_);    // 2/π
        [[nodiscard]] static Float inv_sqrt_pi(int precision = default_precision_);   // 1/√π
        [[nodiscard]] static Float two_inv_sqrt_pi(int precision = default_precision_); // 2/√π

        // --- 平方根・累乗根 ---
        [[nodiscard]] static Float sqrt3(int precision = default_precision_);         // √3
        [[nodiscard]] static Float sqrt5(int precision = default_precision_);         // √5
        [[nodiscard]] static Float inv_sqrt2(int precision = default_precision_);     // 1/√2 = √2/2
        [[nodiscard]] static Float cbrt2(int precision = default_precision_);         // ∛2

        // --- 対数 ---
        [[nodiscard]] static Float ln3(int precision = default_precision_);           // log(3)
        [[nodiscard]] static Float ln5(int precision = default_precision_);           // log(5)
        [[nodiscard]] static Float log2e(int precision = default_precision_);         // log₂(e) = 1/ln2
        [[nodiscard]] static Float log10e(int precision = default_precision_);        // log₁₀(e) = 1/ln10

        // --- 黄金比・基本定数 ---
        [[nodiscard]] static Float phi(int precision = default_precision_);           // φ = (1+√5)/2
        [[nodiscard]] static Float sin1(int precision = default_precision_);          // sin(1)
        [[nodiscard]] static Float cos1(int precision = default_precision_);          // cos(1)
        [[nodiscard]] static Float degree(int precision = default_precision_);        // π/180
        [[nodiscard]] static Float egamma_exp(int precision = default_precision_);    // e^γ

        // --- ゼータ・級数定数 ---
        [[nodiscard]] static Float zeta7(int precision = default_precision_);         // ζ(7)

        // --- 特殊定数 ---
        [[nodiscard]] static Float glaisher(int precision = default_precision_);      // Glaisher-Kinkelin A
        [[nodiscard]] static Float khinchin(int precision = default_precision_);      // Khinchin K
        [[nodiscard]] static Float omega(int precision = default_precision_);         // Ω (Lambert W(1))
        [[nodiscard]] static Float plastic(int precision = default_precision_);       // plastic number ρ
        [[nodiscard]] static Float twin_prime(int precision = default_precision_);    // twin prime C₂
        [[nodiscard]] static Float landau_ramanujan(int precision = default_precision_); // Landau-Ramanujan K
        [[nodiscard]] static Float meissel_mertens(int precision = default_precision_); // Meissel-Mertens M
        [[nodiscard]] static Float bernstein(int precision = default_precision_);     // Bernstein β
        [[nodiscard]] static Float gauss_kuzmin(int precision = default_precision_);  // Gauss-Kuzmin-Wirsing λ
        [[nodiscard]] static Float feigenbaum_delta(int precision = default_precision_); // Feigenbaum δ
        [[nodiscard]] static Float feigenbaum_alpha(int precision = default_precision_); // Feigenbaum α
        [[nodiscard]] static Float erdos_borwein(int precision = default_precision_); // Erdős-Borwein E
        [[nodiscard]] static Float laplace_limit(int precision = default_precision_); // Laplace limit λ*
        [[nodiscard]] static Float soldner(int precision = default_precision_);       // Ramanujan-Soldner μ
        [[nodiscard]] static Float backhouse(int precision = default_precision_);     // Backhouse B
        [[nodiscard]] static Float porter(int precision = default_precision_);        // Porter C
        [[nodiscard]] static Float lieb_square_ice(int precision = default_precision_); // Lieb (8√3/9)
        [[nodiscard]] static Float niven(int precision = default_precision_);         // Niven C
        [[nodiscard]] static Float reciprocal_fibonacci(int precision = default_precision_); // ψ = Σ 1/F(k)
        [[nodiscard]] static Float sierpinski(int precision = default_precision_);    // Sierpiński K
        [[nodiscard]] static Float mills(int precision = default_precision_);         // Mills θ
        [[nodiscard]] static Float dottie(int precision = default_precision_);        // Dottie d
        [[nodiscard]] static Float golomb_dickman(int precision = default_precision_); // Golomb-Dickman λ
        [[nodiscard]] static Float salem(int precision = default_precision_);         // smallest Salem (Lehmer) τ
        [[nodiscard]] static Float cahen(int precision = default_precision_);         // Cahen C
        [[nodiscard]] static Float levy(int precision = default_precision_);          // Lévy β
        [[nodiscard]] static Float copeland_erdos(int precision = default_precision_); // Copeland-Erdős C_CE
        [[nodiscard]] static Float pi_squared_over_6(int precision = default_precision_); // π²/6 = ζ(2)
        [[nodiscard]] static Float pi_squared_over_12(int precision = default_precision_); // π²/12 = ζ(2)/2

        /**
         * @brief Int型に変換
         * @return Int型の値（小数部は切り捨て）
         */
        [[nodiscard]] Int toInt() const;

        /**
         * @brief ビット長を取得
         * @return ビット長
         */
        [[nodiscard]] int bitLength() const;

        /**
         * @brief xを正確に表現するための最小精度（ビット数）を取得
         * @return 最小精度ビット数（0, NaN, Infの場合は0）
         */
        [[nodiscard]] int minPrec() const;

        /**
         * @brief 正の値かどうかを判定
         * @return 正（0以上）ならtrue
         */
        [[nodiscard]] bool isPositive() const { return !isNegative(); }

        /**
         * @brief 整数値かどうかを判定
         * @return 整数ならtrue（NaN/∞はfalse）
         */
        [[nodiscard]] bool isInteger() const;

        /**
         * @brief int型に収まるかを判定（小数部は切り捨て）
         * @return int範囲内ならtrue
         */
        [[nodiscard]] bool fitsInt() const;

        /**
         * @brief int64_t型に収まるかを判定（小数部は切り捨て）
         * @return int64_t範囲内ならtrue
         */
        [[nodiscard]] bool fitsInt64() const;

        /**
         * @brief double型で近似表現可能かを判定（範囲チェック）
         * @return double範囲内ならtrue
         */
        [[nodiscard]] bool fitsDouble() const;

        /**
         * @brief 左シフト代入演算子
         * @param shift シフト量
         * @return 自身への参照
         */
        Float& operator<<=(int shift);

        /**
         * @brief 右シフト代入演算子
         * @param shift シフト量
         * @return 自身への参照
         */
        Float& operator>>=(int shift);

        /**
         * @brief 左シフト演算子
         * @param shift シフト量
         * @return シフト結果
         */
        friend Float operator<<(const Float& value, int shift);

        /**
         * @brief 右シフト演算子
         * @param shift シフト量
         * @return シフト結果
         */
        friend Float operator>>(const Float& value, int shift);

    private:
        // 内部実装ヘルパー関数

        /**
         * @brief 符号なし加算の実装
         * @param rhs 加算する値
         * @return 結果
         */
        // mpn 演算結果から Float を直接構築するタグ型
        struct PreNormalizedTag {};

        // 前提: mantissa は MSB/LSB 正規化済み (normalize()/trimZeroWords() 不要)
        Float(Int&& mantissa, int64_t exponent, bool is_negative, PreNormalizedTag)
            : mantissa_(std::move(mantissa)), exponent_(exponent), is_negative_(is_negative),
              is_infinity_(false), is_nan_(false),
              effective_bits_(INT_MAX), requested_bits_(INT_MAX) {
            checkExponentBounds();
        }

        Float addUnsigned(const Float& rhs) const &;
        Float addUnsigned(const Float& rhs) &&;

        /**
         * @brief 符号なし減算の実装
         * @param rhs 減算する値
         * @return 結果
         */
        Float subtractUnsigned(const Float& rhs) const &;
        Float subtractUnsigned(const Float& rhs) &&;

        // in-place 版: this の mantissa_ バッファを直接上書きする
        // operator+= / -= から呼ばれる。effective_bits_ / requested_bits_ は呼び出し側で設定。
        // is_negative_ は変更しない (addUnsigned) / subtractUnsigned と同じセマンティクス。
        void addUnsignedInPlace(const Float& rhs);
        void subtractUnsignedInPlace(const Float& rhs);

        /**
         * @brief 符号なし比較
         * @param rhs 比較する値
         * @return -1: 自身が小さい, 0: 等しい, 1: 自身が大きい
         */
        int compareUnsigned(const Float& rhs) const;

        /**
         * @brief 丸め処理の実装
         * @param precision_bits 精度（ビット単位）
         * @param mode 丸めモード
         */
        void round(int precision_bits, RoundingMode mode = rounding_mode_);
    };

    // numeric_traits 特殊化の宣言
    template<>
    struct numeric_traits<Float> {
        using value_type = Float;
        using category = floating_point_tag;

        static constexpr bool is_supported = true;
        static constexpr bool is_complex = false;
        static constexpr bool is_integer = false;
        static constexpr bool is_floating_point = true;

        static Float zero() { return Float::zero(); }
        static Float one() { return Float::one(); }
        static Float epsilon() { return Float::epsilon(); }

        static Float abs(const Float& value);
        static Float conj(const Float& value) { return value; }  // 実数の共役は自分自身
        static Float norm(const Float& value);  // |x|^2
        static bool pivotBetter(const Float& a, const Float& b) {
            return abs(a) > abs(b);
        }
    };

    // 基本演算子
    [[nodiscard]] Float operator+(const Float& lhs, const Float& rhs);
    [[nodiscard]] Float operator+(Float&& lhs, const Float& rhs);
    [[nodiscard]] Float operator+(const Float& lhs, Float&& rhs);
    [[nodiscard]] Float operator+(Float&& lhs, Float&& rhs);
    [[nodiscard]] Float operator-(const Float& lhs, const Float& rhs);
    [[nodiscard]] Float operator-(Float&& lhs, const Float& rhs);
    [[nodiscard]] Float operator-(const Float& lhs, Float&& rhs);
    [[nodiscard]] Float operator-(Float&& lhs, Float&& rhs);
    [[nodiscard]] Float operator-(const Float& value);  // 単項マイナス
    [[nodiscard]] Float operator*(const Float& lhs, const Float& rhs);
    [[nodiscard]] Float operator*(Float&& lhs, const Float& rhs);
    [[nodiscard]] Float operator*(const Float& lhs, Float&& rhs);
    [[nodiscard]] Float operator*(Float&& lhs, Float&& rhs);
    [[nodiscard]] Float operator/(const Float& lhs, const Float& rhs);
    [[nodiscard]] Float operator/(Float&& lhs, const Float& rhs);
    // 単一ワード乗除算 — 実装関数 (Float.cpp)
    [[nodiscard]] Float mulScalarF(const Float& lhs, uint64_t rhs);
    [[nodiscard]] Float mulScalarF(Float&& lhs, uint64_t rhs);
    [[nodiscard]] Float mulScalarF(const Float& lhs, int64_t rhs);
    [[nodiscard]] Float mulScalarF(Float&& lhs, int64_t rhs);
    [[nodiscard]] Float divScalarF(const Float& lhs, uint64_t rhs);
    [[nodiscard]] Float divScalarF(Float&& lhs, uint64_t rhs);
    [[nodiscard]] Float divScalarF(const Float& lhs, int64_t rhs);
    [[nodiscard]] Float divScalarF(Float&& lhs, int64_t rhs);

    // テンプレート演算子: Float * scalar, Float / scalar
    template<typename T, std::enable_if_t<std::is_integral_v<T> && !std::is_same_v<T, bool>, int> = 0>
    [[nodiscard]] Float operator*(const Float& lhs, T rhs) {
        if constexpr (std::is_signed_v<T>) return mulScalarF(lhs, static_cast<int64_t>(rhs));
        else return mulScalarF(lhs, static_cast<uint64_t>(rhs));
    }
    template<typename T, std::enable_if_t<std::is_integral_v<T> && !std::is_same_v<T, bool>, int> = 0>
    [[nodiscard]] Float operator*(Float&& lhs, T rhs) {
        if constexpr (std::is_signed_v<T>) return mulScalarF(std::move(lhs), static_cast<int64_t>(rhs));
        else return mulScalarF(std::move(lhs), static_cast<uint64_t>(rhs));
    }
    template<typename T, std::enable_if_t<std::is_integral_v<T> && !std::is_same_v<T, bool>, int> = 0>
    [[nodiscard]] Float operator*(T lhs, const Float& rhs) {
        if constexpr (std::is_signed_v<T>) return mulScalarF(rhs, static_cast<int64_t>(lhs));
        else return mulScalarF(rhs, static_cast<uint64_t>(lhs));
    }
    template<typename T, std::enable_if_t<std::is_integral_v<T> && !std::is_same_v<T, bool>, int> = 0>
    [[nodiscard]] Float operator/(const Float& lhs, T rhs) {
        if constexpr (std::is_signed_v<T>) return divScalarF(lhs, static_cast<int64_t>(rhs));
        else return divScalarF(lhs, static_cast<uint64_t>(rhs));
    }
    template<typename T, std::enable_if_t<std::is_integral_v<T> && !std::is_same_v<T, bool>, int> = 0>
    [[nodiscard]] Float operator/(Float&& lhs, T rhs) {
        if constexpr (std::is_signed_v<T>) return divScalarF(std::move(lhs), static_cast<int64_t>(rhs));
        else return divScalarF(std::move(lhs), static_cast<uint64_t>(rhs));
    }

    // double との混合演算 — Float(double) は explicit なので演算子を提供
    // float 版は不要: float の仮数部 (23 bit) は 1 リムに収まり、
    // Float(float) と Float(double) は同じ内部表現になる。
    // C++ の暗黙変換 float→double により、ここの double 版が使われる。
    [[nodiscard]] inline Float operator+(const Float& lhs, double rhs) { return lhs + Float(rhs); }
    [[nodiscard]] inline Float operator+(Float&& lhs, double rhs)      { return std::move(lhs) + Float(rhs); }
    [[nodiscard]] inline Float operator+(double lhs, const Float& rhs) { return Float(lhs) + rhs; }
    [[nodiscard]] inline Float operator-(const Float& lhs, double rhs) { return lhs - Float(rhs); }
    [[nodiscard]] inline Float operator-(Float&& lhs, double rhs)      { return std::move(lhs) - Float(rhs); }
    [[nodiscard]] inline Float operator-(double lhs, const Float& rhs) { return Float(lhs) - rhs; }
    [[nodiscard]] inline Float operator*(const Float& lhs, double rhs) { return lhs * Float(rhs); }
    [[nodiscard]] inline Float operator*(Float&& lhs, double rhs)      { return std::move(lhs) * Float(rhs); }
    [[nodiscard]] inline Float operator*(double lhs, const Float& rhs) { return Float(lhs) * rhs; }
    [[nodiscard]] inline Float operator/(const Float& lhs, double rhs) { return lhs / Float(rhs); }
    [[nodiscard]] inline Float operator/(Float&& lhs, double rhs)      { return std::move(lhs) / Float(rhs); }
    [[nodiscard]] inline Float operator/(double lhs, const Float& rhs) { return Float(lhs) / rhs; }
    [[nodiscard]] inline bool operator==(const Float& lhs, double rhs) { return lhs == Float(rhs); }
    [[nodiscard]] inline bool operator==(double lhs, const Float& rhs) { return Float(lhs) == rhs; }
    [[nodiscard]] inline std::partial_ordering operator<=>(const Float& lhs, double rhs) { return lhs <=> Float(rhs); }
    [[nodiscard]] inline std::partial_ordering operator<=>(double lhs, const Float& rhs) { return Float(lhs) <=> rhs; }

    [[nodiscard]] std::partial_ordering operator<=>(const Float& lhs, const Float& rhs);
    [[nodiscard]] bool operator==(const Float& lhs, const Float& rhs);

    // ストリーム入出力
    std::ostream& operator<<(std::ostream& os, const Float& value);
    std::istream& operator>>(std::istream& is, Float& value);

    // 浮動小数点ユーティリティ関数
    [[nodiscard]] Float fmin(const Float& a, const Float& b);
    [[nodiscard]] Float fmin(Float&& a, Float&& b);
    [[nodiscard]] Float fmax(const Float& a, const Float& b);
    [[nodiscard]] Float fmax(Float&& a, Float&& b);
    [[nodiscard]] Float fdim(const Float& a, const Float& b);
    [[nodiscard]] Float fdim(Float&& a, Float&& b);
    [[nodiscard]] Float copySign(const Float& x, const Float& y);
    [[nodiscard]] Float copySign(Float&& x, const Float& y);
    [[nodiscard]] bool signBit(const Float& x);

    // 線形補間・中点
    [[nodiscard]] Float lerp(const Float& a, const Float& b, const Float& t, int precision);
    [[nodiscard]] Float midpoint(const Float& a, const Float& b);

    // <cmath> 互換ユーティリティ
    [[nodiscard]] Float modf(const Float& x, Float& iptr);        // 整数部+小数部分離
    [[nodiscard]] int64_t ilogb(const Float& x);                   // floor(log2(|x|))
    [[nodiscard]] Float logb(const Float& x);                      // floor(log2(|x|)) as Float
    [[nodiscard]] Float scalbn(const Float& x, int n);             // x * 2^n (= ldexp)
    [[nodiscard]] Float scalbn(Float&& x, int n);
    [[nodiscard]] Float nearbyint(const Float& x);                 // round に同じ
    [[nodiscard]] Float nearbyint(Float&& x);
    [[nodiscard]] Float rint(const Float& x);                      // round に同じ
    [[nodiscard]] Float rint(Float&& x);

    // 丸め・整数化関数
    [[nodiscard]] Float floor(const Float& x);
    [[nodiscard]] Float floor(Float&& x);
    [[nodiscard]] Float ceil(const Float& x);
    [[nodiscard]] Float ceil(Float&& x);
    [[nodiscard]] Float trunc(const Float& x);
    [[nodiscard]] Float trunc(Float&& x);
    [[nodiscard]] Float round(const Float& x);
    [[nodiscard]] Float round(Float&& x);
    [[nodiscard]] Float roundEven(const Float& x);               // 最近接偶数丸め (banker's rounding)
    [[nodiscard]] Float roundEven(Float&& x);
    [[nodiscard]] Float frac(const Float& x);
    [[nodiscard]] Float frac(Float&& x);

    // 逆三角関数
    [[nodiscard]] Float asin(const Float& x, int precision);
    [[nodiscard]] Float asin(Float&& x, int precision);
    [[nodiscard]] Float acos(const Float& x, int precision);
    [[nodiscard]] Float acos(Float&& x, int precision);
    [[nodiscard]] Float atan(const Float& x, int precision);
    [[nodiscard]] Float atan(Float&& x, int precision);
    [[nodiscard]] Float atan2(const Float& y, const Float& x, int precision);
    [[nodiscard]] Float atan2(Float&& y, Float&& x, int precision);

    // 逆三角関数 (π 単位の結果)
    [[nodiscard]] Float asinPi(const Float& x, int precision);
    [[nodiscard]] Float asinPi(Float&& x, int precision);
    [[nodiscard]] Float acosPi(const Float& x, int precision);
    [[nodiscard]] Float acosPi(Float&& x, int precision);
    [[nodiscard]] Float atanPi(const Float& x, int precision);
    [[nodiscard]] Float atanPi(Float&& x, int precision);
    [[nodiscard]] Float atan2Pi(const Float& y, const Float& x, int precision);
    [[nodiscard]] Float atan2Pi(Float&& y, Float&& x, int precision);

    // 逆双曲線関数
    [[nodiscard]] Float asinh(const Float& x, int precision);
    [[nodiscard]] Float asinh(Float&& x, int precision);
    [[nodiscard]] Float acosh(const Float& x, int precision);
    [[nodiscard]] Float acosh(Float&& x, int precision);
    [[nodiscard]] Float atanh(const Float& x, int precision);
    [[nodiscard]] Float atanh(Float&& x, int precision);

    // 対数・指数バリエーション
    [[nodiscard]] Float log2(const Float& x, int precision);
    [[nodiscard]] Float log2(Float&& x, int precision);
    [[nodiscard]] Float log10(const Float& x, int precision);
    [[nodiscard]] Float log10(Float&& x, int precision);
    [[nodiscard]] Float log1p(const Float& x, int precision);
    [[nodiscard]] Float log1p(Float&& x, int precision);
    [[nodiscard]] Float exp2(const Float& x, int precision);
    [[nodiscard]] Float exp2(Float&& x, int precision);
    [[nodiscard]] Float exp10(const Float& x, int precision);
    [[nodiscard]] Float exp10(Float&& x, int precision);
    [[nodiscard]] Float expm1(const Float& x, int precision);
    [[nodiscard]] Float expm1(Float&& x, int precision);
    [[nodiscard]] Float exp2m1(const Float& x, int precision);
    [[nodiscard]] Float exp2m1(Float&& x, int precision);
    [[nodiscard]] Float exp10m1(const Float& x, int precision);
    [[nodiscard]] Float exp10m1(Float&& x, int precision);
    [[nodiscard]] Float log2p1(const Float& x, int precision);
    [[nodiscard]] Float log2p1(Float&& x, int precision);
    [[nodiscard]] Float log10p1(const Float& x, int precision);
    [[nodiscard]] Float log10p1(Float&& x, int precision);
    [[nodiscard]] Float compound(const Float& x, int n, int precision);
    [[nodiscard]] Float compound(Float&& x, int n, int precision);

    // n乗根・逆平方根
    [[nodiscard]] Float cbrt(const Float& x, int precision);
    [[nodiscard]] Float cbrt(Float&& x, int precision);
    [[nodiscard]] Float nthRoot(const Float& x, int n, int precision);
    [[nodiscard]] Float nthRoot(Float&& x, int n, int precision);
    [[nodiscard]] Float recSqrt(const Float& x, int precision);
    [[nodiscard]] Float recSqrt(Float&& x, int precision);

    // 融合乗加算
    [[nodiscard]] Float fma(const Float& a, const Float& b, const Float& c, int precision);
    [[nodiscard]] Float fma(Float&& a, Float&& b, Float&& c, int precision);
    [[nodiscard]] Float fms(const Float& a, const Float& b, const Float& c, int precision);
    [[nodiscard]] Float fms(Float&& a, Float&& b, Float&& c, int precision);

    // 二重積和/差 (MPFR mpfr_fmma/mpfr_fmms 相当)
    [[nodiscard]] Float fmma(const Float& a, const Float& b, const Float& c, const Float& d, int precision);
    [[nodiscard]] Float fmms(const Float& a, const Float& b, const Float& c, const Float& d, int precision);

    // fmod / remainder / hypot
    [[nodiscard]] Float fmod(const Float& x, const Float& y);
    [[nodiscard]] Float fmod(Float&& x, Float&& y);
    [[nodiscard]] Float remainder(const Float& x, const Float& y);
    [[nodiscard]] Float remainder(Float&& x, Float&& y);
    [[nodiscard]] std::pair<Float, int> remquo(const Float& x, const Float& y);
    [[nodiscard]] std::pair<Float, int> remquo(Float&& x, Float&& y);
    [[nodiscard]] Float hypot(const Float& x, const Float& y, int precision);
    [[nodiscard]] Float hypot(Float&& x, Float&& y, int precision);

    // 同時計算
    void sinCos(const Float& x, Float& sin_result, Float& cos_result, int precision);
    void sinCos(Float&& x, Float& sin_result, Float& cos_result, int precision);
    void sinhCosh(const Float& x, Float& sinh_result, Float& cosh_result, int precision);
    void sinhCosh(Float&& x, Float& sinh_result, Float& cosh_result, int precision);

    // 逆数三角関数・双曲線関数
    [[nodiscard]] Float sec(const Float& x, int precision);
    [[nodiscard]] Float sec(Float&& x, int precision);
    [[nodiscard]] Float csc(const Float& x, int precision);
    [[nodiscard]] Float csc(Float&& x, int precision);
    [[nodiscard]] Float cot(const Float& x, int precision);
    [[nodiscard]] Float cot(Float&& x, int precision);
    [[nodiscard]] Float sech(const Float& x, int precision);
    [[nodiscard]] Float sech(Float&& x, int precision);
    [[nodiscard]] Float csch(const Float& x, int precision);
    [[nodiscard]] Float csch(Float&& x, int precision);
    [[nodiscard]] Float coth(const Float& x, int precision);
    [[nodiscard]] Float coth(Float&& x, int precision);

    // 階乗
    [[nodiscard]] Float factorial(int n, int precision);

    // πベース三角関数
    [[nodiscard]] Float sinPi(const Float& x, int precision);
    [[nodiscard]] Float sinPi(Float&& x, int precision);
    [[nodiscard]] Float cosPi(const Float& x, int precision);
    [[nodiscard]] Float cosPi(Float&& x, int precision);
    [[nodiscard]] Float tanPi(const Float& x, int precision);
    [[nodiscard]] Float tanPi(Float&& x, int precision);

    // 任意角度単位三角関数 (IEEE 754-2019)
    // sinu(x, u) = sin(2πx/u), cosu(x, u) = cos(2πx/u), tanu(x, u) = tan(2πx/u)
    // u は全回転の単位数 (360=度, 400=グラジアン)
    [[nodiscard]] Float sinu(const Float& x, int u, int precision);
    [[nodiscard]] Float sinu(Float&& x, int u, int precision);
    [[nodiscard]] Float cosu(const Float& x, int u, int precision);
    [[nodiscard]] Float cosu(Float&& x, int u, int precision);
    [[nodiscard]] Float tanu(const Float& x, int u, int precision);
    [[nodiscard]] Float tanu(Float&& x, int u, int precision);

    // nextAbove / nextBelow
    [[nodiscard]] Float nextAbove(const Float& x);
    [[nodiscard]] Float nextAbove(Float&& x);
    [[nodiscard]] Float nextBelow(const Float& x);
    [[nodiscard]] Float nextBelow(Float&& x);

    // 丸め方向判定 (MPFR mpfr_can_round 互換)
    // x: err_bits ビットの精度で rnd1 丸めにより計算された近似値
    // target_prec: 目標精度 (10進桁数)
    // rnd2 で target_prec に丸めたとき、丸め方向が一意に確定するなら true
    [[nodiscard]] bool canRound(const Float& x, int err_bits,
                                RoundingMode rnd1, RoundingMode rnd2,
                                int target_prec);

    // AGM (算術幾何平均)
    [[nodiscard]] Float agm(const Float& a, const Float& b, int precision);
    [[nodiscard]] Float agm(Float&& a, Float&& b, int precision);

    // 高精度総和・内積
    [[nodiscard]] Float sum(std::span<const Float> values, int precision);
    [[nodiscard]] Float dot(std::span<const Float> a, std::span<const Float> b, int precision);

    // 誤差関数
    [[nodiscard]] Float erf(const Float& x, int precision);
    [[nodiscard]] Float erf(Float&& x, int precision);
    [[nodiscard]] Float erfc(const Float& x, int precision);
    [[nodiscard]] Float erfc(Float&& x, int precision);
    [[nodiscard]] Float erfcx(const Float& x, int precision);    // exp(x²)·erfc(x)
    [[nodiscard]] Float erfcx(Float&& x, int precision);

    // ガンマ関数・関連関数
    [[nodiscard]] Float lnGamma(const Float& x, int precision);
    [[nodiscard]] Float lnGamma(Float&& x, int precision);
    [[nodiscard]] Float gamma(const Float& x, int precision);
    [[nodiscard]] Float gamma(Float&& x, int precision);
    [[nodiscard]] Float beta(const Float& a, const Float& b, int precision);
    [[nodiscard]] Float beta(Float&& a, Float&& b, int precision);
    [[nodiscard]] Float digamma(const Float& x, int precision);
    [[nodiscard]] Float digamma(Float&& x, int precision);
    [[nodiscard]] Float trigamma(const Float& x, int precision);
    [[nodiscard]] Float trigamma(Float&& x, int precision);
    [[nodiscard]] Float polygamma(int n, const Float& x, int precision);
    [[nodiscard]] Float polygamma(int n, Float&& x, int precision);

    // 不完全ガンマ関数
    [[nodiscard]] Float gammaP(const Float& a, const Float& x, int precision);     // P(a,x) = γ(a,x)/Γ(a) 正規化下側
    [[nodiscard]] Float gammaP(Float&& a, Float&& x, int precision);
    [[nodiscard]] Float gammaQ(const Float& a, const Float& x, int precision);     // Q(a,x) = Γ(a,x)/Γ(a) 正規化上側
    [[nodiscard]] Float gammaQ(Float&& a, Float&& x, int precision);
    [[nodiscard]] Float gammaLower(const Float& a, const Float& x, int precision); // γ(a,x) 下側不完全ガンマ
    [[nodiscard]] Float gammaLower(Float&& a, Float&& x, int precision);
    [[nodiscard]] Float gammaUpper(const Float& a, const Float& x, int precision); // Γ(a,x) 上側不完全ガンマ
    [[nodiscard]] Float gammaUpper(Float&& a, Float&& x, int precision);

    // 不完全ベータ関数
    [[nodiscard]] Float betaRegularized(const Float& x, const Float& a, const Float& b, int precision); // I_x(a,b)
    [[nodiscard]] Float betaRegularized(Float&& x, Float&& a, Float&& b, int precision);

    // 楕円積分 — Carlson 対称形式
    [[nodiscard]] Float carlsonRF(const Float& x, const Float& y, const Float& z, int precision);
    [[nodiscard]] Float carlsonRF(Float&& x, Float&& y, Float&& z, int precision);
    [[nodiscard]] Float carlsonRC(const Float& x, const Float& y, int precision);
    [[nodiscard]] Float carlsonRC(Float&& x, Float&& y, int precision);
    [[nodiscard]] Float carlsonRD(const Float& x, const Float& y, const Float& z, int precision);
    [[nodiscard]] Float carlsonRD(Float&& x, Float&& y, Float&& z, int precision);
    [[nodiscard]] Float carlsonRJ(const Float& x, const Float& y, const Float& z, const Float& p, int precision);
    [[nodiscard]] Float carlsonRJ(Float&& x, Float&& y, Float&& z, Float&& p, int precision);

    // 楕円積分 — Legendre 形式 (完全)
    [[nodiscard]] Float ellipticK(const Float& k, int precision);                                      // K(k) 第1種
    [[nodiscard]] Float ellipticK(Float&& k, int precision);
    [[nodiscard]] Float ellipticE(const Float& k, int precision);                                      // E(k) 第2種
    [[nodiscard]] Float ellipticE(Float&& k, int precision);
    [[nodiscard]] Float ellipticPi(const Float& n, const Float& k, int precision);                     // Π(n,k) 第3種
    [[nodiscard]] Float ellipticPi(Float&& n, Float&& k, int precision);

    // 楕円積分 — Legendre 形式 (不完全)
    [[nodiscard]] Float ellipticF(const Float& phi, const Float& k, int precision);                    // F(φ,k) 第1種
    [[nodiscard]] Float ellipticF(Float&& phi, Float&& k, int precision);
    [[nodiscard]] Float ellipticE(const Float& phi, const Float& k, int precision);                    // E(φ,k) 第2種
    [[nodiscard]] Float ellipticE(Float&& phi, Float&& k, int precision);
    [[nodiscard]] Float ellipticPi(const Float& n, const Float& phi, const Float& k, int precision);   // Π(n,φ,k) 第3種
    [[nodiscard]] Float ellipticPi(Float&& n, Float&& phi, Float&& k, int precision);

    // Jacobi 楕円関数
    [[nodiscard]] Float jacobiSn(const Float& u, const Float& k, int precision);
    [[nodiscard]] Float jacobiSn(Float&& u, Float&& k, int precision);
    [[nodiscard]] Float jacobiCn(const Float& u, const Float& k, int precision);
    [[nodiscard]] Float jacobiCn(Float&& u, Float&& k, int precision);
    [[nodiscard]] Float jacobiDn(const Float& u, const Float& k, int precision);
    [[nodiscard]] Float jacobiDn(Float&& u, Float&& k, int precision);

    // ゼータ関数
    [[nodiscard]] Float zeta(const Float& s, int precision);
    [[nodiscard]] Float zeta(Float&& s, int precision);
    [[nodiscard]] Float hurwitzZeta(const Float& s, const Float& a, int precision);  // ζ(s,a)
    [[nodiscard]] Float hurwitzZeta(Float&& s, Float&& a, int precision);
    [[nodiscard]] Float dirichletEta(const Float& s, int precision);                 // η(s)
    [[nodiscard]] Float dirichletEta(Float&& s, int precision);

    // 指数積分・ディロガリズム
    [[nodiscard]] Float expint(const Float& x, int precision);   // Ei(x)
    [[nodiscard]] Float expint(Float&& x, int precision);
    [[nodiscard]] Float li(const Float& x, int precision);        // li(x) = Ei(ln(x))
    [[nodiscard]] Float li(Float&& x, int precision);
    [[nodiscard]] Float dilog(const Float& x, int precision);     // Li₂(x)
    [[nodiscard]] Float dilog(Float&& x, int precision);

    // ベッセル関数
    [[nodiscard]] Float besselJ(int n, const Float& x, int precision);  // J_n(x)
    [[nodiscard]] Float besselJ(int n, Float&& x, int precision);
    [[nodiscard]] Float besselY(int n, const Float& x, int precision);  // Y_n(x)
    [[nodiscard]] Float besselY(int n, Float&& x, int precision);
    [[nodiscard]] Float besselI(int n, const Float& x, int precision);  // I_n(x) 変形第1種
    [[nodiscard]] Float besselI(int n, Float&& x, int precision);
    [[nodiscard]] Float besselK(int n, const Float& x, int precision);  // K_n(x) 変形第2種
    [[nodiscard]] Float besselK(int n, Float&& x, int precision);
    [[nodiscard]] Float sphericalBesselJ(int n, const Float& x, int precision);  // j_n(x) 球ベッセル第1種
    [[nodiscard]] Float sphericalBesselJ(int n, Float&& x, int precision);
    [[nodiscard]] Float sphericalBesselY(int n, const Float& x, int precision);  // y_n(x) 球ベッセル第2種
    [[nodiscard]] Float sphericalBesselY(int n, Float&& x, int precision);

    // エアリー関数
    [[nodiscard]] Float airyAi(const Float& x, int precision);       // Ai(x)
    [[nodiscard]] Float airyAi(Float&& x, int precision);
    [[nodiscard]] Float airyBi(const Float& x, int precision);       // Bi(x)
    [[nodiscard]] Float airyBi(Float&& x, int precision);
    [[nodiscard]] Float airyAiPrime(const Float& x, int precision);  // Ai'(x)
    [[nodiscard]] Float airyAiPrime(Float&& x, int precision);
    [[nodiscard]] Float airyBiPrime(const Float& x, int precision);  // Bi'(x)
    [[nodiscard]] Float airyBiPrime(Float&& x, int precision);

    // 超幾何関数
    [[nodiscard]] Float confHyperg(const Float& a, const Float& b, const Float& z, int precision);                         // ₁F₁(a; b; z)
    [[nodiscard]] Float confHyperg(Float&& a, Float&& b, Float&& z, int precision);
    [[nodiscard]] Float hyperg(const Float& a, const Float& b, const Float& c, const Float& z, int precision);             // ₂F₁(a, b; c; z)
    [[nodiscard]] Float hyperg(Float&& a, Float&& b, Float&& c, Float&& z, int precision);
    [[nodiscard]] Float hyperg0F1(const Float& b, const Float& z, int precision);                                          // ₀F₁(; b; z)
    [[nodiscard]] Float hyperg0F1(Float&& b, Float&& z, int precision);

    // Legendre 多項式
    [[nodiscard]] Float legendreP(int n, const Float& x, int precision);              // P_n(x)
    [[nodiscard]] Float legendreP(int n, Float&& x, int precision);
    [[nodiscard]] Float assocLegendreP(int n, int m, const Float& x, int precision);  // P_n^m(x)
    [[nodiscard]] Float assocLegendreP(int n, int m, Float&& x, int precision);

    // Hermite 多項式 (物理学者版: H_n(x))
    [[nodiscard]] Float hermite(int n, const Float& x, int precision);
    [[nodiscard]] Float hermite(int n, Float&& x, int precision);

    // Laguerre 多項式
    [[nodiscard]] Float laguerre(int n, const Float& x, int precision);               // L_n(x)
    [[nodiscard]] Float laguerre(int n, Float&& x, int precision);
    [[nodiscard]] Float assocLaguerre(int n, int m, const Float& x, int precision);   // L_n^m(x)
    [[nodiscard]] Float assocLaguerre(int n, int m, Float&& x, int precision);

    // Lambert W 関数
    [[nodiscard]] Float lambertW0(const Float& x, int precision);    // W₀(x)
    [[nodiscard]] Float lambertW0(Float&& x, int precision);
    [[nodiscard]] Float lambertWm1(const Float& x, int precision);   // W₋₁(x)
    [[nodiscard]] Float lambertWm1(Float&& x, int precision);

    // 指数積分・三角積分
    [[nodiscard]] Float expintN(int n, const Float& x, int precision);        // E_n(x)
    [[nodiscard]] Float expintN(int n, Float&& x, int precision);
    [[nodiscard]] Float sinIntegral(const Float& x, int precision);           // Si(x)
    [[nodiscard]] Float sinIntegral(Float&& x, int precision);
    [[nodiscard]] Float cosIntegral(const Float& x, int precision);           // Ci(x)
    [[nodiscard]] Float cosIntegral(Float&& x, int precision);

    // 乱数生成
    [[nodiscard]] Float randomFloat(int precision);                                    // [0, 1) 一様分布
    [[nodiscard]] Float randomFloat(const Float& lo, const Float& hi, int precision);  // [lo, hi) 一様分布
    [[nodiscard]] Float randomFloat(Float&& lo, Float&& hi, int precision);
    [[nodiscard]] Float normalRandom(int precision);                                   // N(0,1) 正規分布
    [[nodiscard]] Float exponentialRandom(int precision);                              // Exp(1) 指数分布

    // ================================================================
    // バイナリシリアライズ
    // フォーマット: [tag:1B][sign:1B][exponent:8B LE][eff:4B LE][req:4B LE][mantissa...]
    //   tag: 0x00=zero, 0x01=normal, 0x02=infinity, 0x03=nan
    //   mantissa: Int::exportBinary の出力 (tag==0x01 のみ)
    // ================================================================
    [[nodiscard]] std::vector<uint8_t> exportBinary(const Float& value);
    [[nodiscard]] Float importBinaryFloat(std::span<const uint8_t> data);

    // ================================================================
    // 精度省略オーバーロード
    // 入力の requestedBits() から精度を決定。exact (INT_MAX) の場合は
    // defaultPrecision() にフォールバック。
    // Complex<Float> の内部から ADL で呼ばれることを想定。
    // ================================================================

    [[nodiscard]] inline Float exp(const Float& x)  { return exp(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float exp(Float&& x)       { int p = Float::requestedPrecision(x); return exp(std::move(x), p); }
    [[nodiscard]] inline Float log(const Float& x)  { return log(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float log(Float&& x)       { int p = Float::requestedPrecision(x); return log(std::move(x), p); }
    [[nodiscard]] inline Float logUi(unsigned long long n) { return logUi(n, Float::defaultPrecision()); }
    [[nodiscard]] inline Float sin(const Float& x)  { return sin(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float sin(Float&& x)       { int p = Float::requestedPrecision(x); return sin(std::move(x), p); }
    [[nodiscard]] inline Float cos(const Float& x)  { return cos(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float cos(Float&& x)       { int p = Float::requestedPrecision(x); return cos(std::move(x), p); }
    [[nodiscard]] inline Float tan(const Float& x)  { return tan(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float tan(Float&& x)       { int p = Float::requestedPrecision(x); return tan(std::move(x), p); }
    [[nodiscard]] inline Float sqr(const Float& x)   { return sqr(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float sqr(Float&& x)       { int p = Float::requestedPrecision(x); return sqr(std::move(x), p); }
    [[nodiscard]] inline Float sqrt(const Float& x) { return sqrt(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float sqrt(Float&& x)      { int p = Float::requestedPrecision(x); return sqrt(std::move(x), p); }
    [[nodiscard]] inline Float sinh(const Float& x) { return sinh(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float sinh(Float&& x)      { int p = Float::requestedPrecision(x); return sinh(std::move(x), p); }
    [[nodiscard]] inline Float cosh(const Float& x) { return cosh(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float cosh(Float&& x)      { int p = Float::requestedPrecision(x); return cosh(std::move(x), p); }
    [[nodiscard]] inline Float tanh(const Float& x) { return tanh(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float tanh(Float&& x)      { int p = Float::requestedPrecision(x); return tanh(std::move(x), p); }
    [[nodiscard]] inline Float asin(const Float& x) { return asin(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float asin(Float&& x)      { int p = Float::requestedPrecision(x); return asin(std::move(x), p); }
    [[nodiscard]] inline Float acos(const Float& x) { return acos(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float acos(Float&& x)      { int p = Float::requestedPrecision(x); return acos(std::move(x), p); }
    [[nodiscard]] inline Float atan(const Float& x) { return atan(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float atan(Float&& x)      { int p = Float::requestedPrecision(x); return atan(std::move(x), p); }
    [[nodiscard]] inline Float atan2(const Float& y, const Float& x) {
        int ry = y.requestedBits(), rx = x.requestedBits();
        int req = (ry >= INT_MAX && rx >= INT_MAX) ? INT_MAX
                : (ry >= INT_MAX) ? rx : (rx >= INT_MAX) ? ry : std::max(ry, rx);
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return atan2(y, x, prec);
    }
    [[nodiscard]] inline Float atan2(Float&& y, Float&& x) {
        int ry = y.requestedBits(), rx = x.requestedBits();
        int req = (ry >= INT_MAX && rx >= INT_MAX) ? INT_MAX
                : (ry >= INT_MAX) ? rx : (rx >= INT_MAX) ? ry : std::max(ry, rx);
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return atan2(std::move(y), std::move(x), prec);
    }
    [[nodiscard]] inline Float asinPi(const Float& x) { return asinPi(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float asinPi(Float&& x)     { int p = Float::requestedPrecision(x); return asinPi(std::move(x), p); }
    [[nodiscard]] inline Float acosPi(const Float& x) { return acosPi(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float acosPi(Float&& x)     { int p = Float::requestedPrecision(x); return acosPi(std::move(x), p); }
    [[nodiscard]] inline Float atanPi(const Float& x) { return atanPi(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float atanPi(Float&& x)     { int p = Float::requestedPrecision(x); return atanPi(std::move(x), p); }
    [[nodiscard]] inline Float atan2Pi(const Float& y, const Float& x) {
        int ry = y.requestedBits(), rx = x.requestedBits();
        int req = (ry >= INT_MAX && rx >= INT_MAX) ? INT_MAX
                : (ry >= INT_MAX) ? rx : (rx >= INT_MAX) ? ry : std::max(ry, rx);
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return atan2Pi(y, x, prec);
    }
    [[nodiscard]] inline Float atan2Pi(Float&& y, Float&& x) {
        int ry = y.requestedBits(), rx = x.requestedBits();
        int req = (ry >= INT_MAX && rx >= INT_MAX) ? INT_MAX
                : (ry >= INT_MAX) ? rx : (rx >= INT_MAX) ? ry : std::max(ry, rx);
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return atan2Pi(std::move(y), std::move(x), prec);
    }
    [[nodiscard]] inline Float asinh(const Float& x) { return asinh(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float asinh(Float&& x)      { int p = Float::requestedPrecision(x); return asinh(std::move(x), p); }
    [[nodiscard]] inline Float acosh(const Float& x) { return acosh(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float acosh(Float&& x)      { int p = Float::requestedPrecision(x); return acosh(std::move(x), p); }
    [[nodiscard]] inline Float atanh(const Float& x) { return atanh(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float atanh(Float&& x)      { int p = Float::requestedPrecision(x); return atanh(std::move(x), p); }
    [[nodiscard]] inline Float pow(const Float& x, const Float& y) {
        int ry = y.requestedBits(), rx = x.requestedBits();
        int req = (ry >= INT_MAX && rx >= INT_MAX) ? INT_MAX
                : (ry >= INT_MAX) ? rx : (rx >= INT_MAX) ? ry : std::max(ry, rx);
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return pow(x, y, prec);
    }
    [[nodiscard]] inline Float pow(Float&& x, Float&& y) {
        int ry = y.requestedBits(), rx = x.requestedBits();
        int req = (ry >= INT_MAX && rx >= INT_MAX) ? INT_MAX
                : (ry >= INT_MAX) ? rx : (rx >= INT_MAX) ? ry : std::max(ry, rx);
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return pow(std::move(x), std::move(y), prec);
    }
    [[nodiscard]] inline Float pow(const Float& x, int n) {
        return pow(x, n, Float::requestedPrecision(x));
    }
    [[nodiscard]] inline Float pow(Float&& x, int n) {
        int p = Float::requestedPrecision(x); return pow(std::move(x), n, p);
    }
    [[nodiscard]] inline Float log2(const Float& x)  { return log2(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float log2(Float&& x)       { int p = Float::requestedPrecision(x); return log2(std::move(x), p); }
    [[nodiscard]] inline Float log10(const Float& x) { return log10(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float log10(Float&& x)      { int p = Float::requestedPrecision(x); return log10(std::move(x), p); }
    [[nodiscard]] inline Float log1p(const Float& x) { return log1p(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float log1p(Float&& x)      { int p = Float::requestedPrecision(x); return log1p(std::move(x), p); }
    [[nodiscard]] inline Float exp2(const Float& x)  { return exp2(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float exp2(Float&& x)       { int p = Float::requestedPrecision(x); return exp2(std::move(x), p); }
    [[nodiscard]] inline Float exp10(const Float& x) { return exp10(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float exp10(Float&& x)      { int p = Float::requestedPrecision(x); return exp10(std::move(x), p); }
    [[nodiscard]] inline Float expm1(const Float& x)  { return expm1(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float expm1(Float&& x)       { int p = Float::requestedPrecision(x); return expm1(std::move(x), p); }
    [[nodiscard]] inline Float exp2m1(const Float& x)  { return exp2m1(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float exp2m1(Float&& x)       { int p = Float::requestedPrecision(x); return exp2m1(std::move(x), p); }
    [[nodiscard]] inline Float exp10m1(const Float& x) { return exp10m1(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float exp10m1(Float&& x)      { int p = Float::requestedPrecision(x); return exp10m1(std::move(x), p); }
    [[nodiscard]] inline Float log2p1(const Float& x)  { return log2p1(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float log2p1(Float&& x)       { int p = Float::requestedPrecision(x); return log2p1(std::move(x), p); }
    [[nodiscard]] inline Float log10p1(const Float& x) { return log10p1(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float log10p1(Float&& x)      { int p = Float::requestedPrecision(x); return log10p1(std::move(x), p); }
    [[nodiscard]] inline Float compound(const Float& x, int n) { return compound(x, n, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float compound(Float&& x, int n)      { int p = Float::requestedPrecision(x); return compound(std::move(x), n, p); }
    [[nodiscard]] inline Float cbrt(const Float& x)  { return cbrt(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float cbrt(Float&& x)       { int p = Float::requestedPrecision(x); return cbrt(std::move(x), p); }
    [[nodiscard]] inline Float nthRoot(const Float& x, int n) { return nthRoot(x, n, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float nthRoot(Float&& x, int n)      { int p = Float::requestedPrecision(x); return nthRoot(std::move(x), n, p); }
    [[nodiscard]] inline Float recSqrt(const Float& x) { return recSqrt(x, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float recSqrt(Float&& x)      { int p = Float::requestedPrecision(x); return recSqrt(std::move(x), p); }
    [[nodiscard]] inline Float hypot(const Float& x, const Float& y) {
        int ry = y.requestedBits(), rx = x.requestedBits();
        int req = (ry >= INT_MAX && rx >= INT_MAX) ? INT_MAX
                : (ry >= INT_MAX) ? rx : (rx >= INT_MAX) ? ry : std::max(ry, rx);
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return hypot(x, y, prec);
    }
    [[nodiscard]] inline Float hypot(Float&& x, Float&& y) {
        int ry = y.requestedBits(), rx = x.requestedBits();
        int req = (ry >= INT_MAX && rx >= INT_MAX) ? INT_MAX
                : (ry >= INT_MAX) ? rx : (rx >= INT_MAX) ? ry : std::max(ry, rx);
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return hypot(std::move(x), std::move(y), prec);
    }
    [[nodiscard]] inline Float fma(const Float& a, const Float& b, const Float& c) {
        int ra = a.requestedBits(), rb = b.requestedBits(), rc = c.requestedBits();
        int req = INT_MAX;
        if (ra < INT_MAX) req = ra;
        if (rb < INT_MAX) req = (req < INT_MAX) ? std::max(req, rb) : rb;
        if (rc < INT_MAX) req = (req < INT_MAX) ? std::max(req, rc) : rc;
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return fma(a, b, c, prec);
    }
    [[nodiscard]] inline Float fma(Float&& a, Float&& b, Float&& c) {
        int ra = a.requestedBits(), rb = b.requestedBits(), rc = c.requestedBits();
        int req = INT_MAX;
        if (ra < INT_MAX) req = ra;
        if (rb < INT_MAX) req = (req < INT_MAX) ? std::max(req, rb) : rb;
        if (rc < INT_MAX) req = (req < INT_MAX) ? std::max(req, rc) : rc;
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return fma(std::move(a), std::move(b), std::move(c), prec);
    }
    [[nodiscard]] inline Float fms(const Float& a, const Float& b, const Float& c) {
        int ra = a.requestedBits(), rb = b.requestedBits(), rc = c.requestedBits();
        int req = INT_MAX;
        if (ra < INT_MAX) req = ra;
        if (rb < INT_MAX) req = (req < INT_MAX) ? std::max(req, rb) : rb;
        if (rc < INT_MAX) req = (req < INT_MAX) ? std::max(req, rc) : rc;
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return fms(a, b, c, prec);
    }
    [[nodiscard]] inline Float fms(Float&& a, Float&& b, Float&& c) {
        int ra = a.requestedBits(), rb = b.requestedBits(), rc = c.requestedBits();
        int req = INT_MAX;
        if (ra < INT_MAX) req = ra;
        if (rb < INT_MAX) req = (req < INT_MAX) ? std::max(req, rb) : rb;
        if (rc < INT_MAX) req = (req < INT_MAX) ? std::max(req, rc) : rc;
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return fms(std::move(a), std::move(b), std::move(c), prec);
    }
    [[nodiscard]] inline Float fmma(const Float& a, const Float& b, const Float& c, const Float& d) {
        int ra = a.requestedBits(), rb = b.requestedBits();
        int rc = c.requestedBits(), rd = d.requestedBits();
        int req = INT_MAX;
        for (int v : {ra, rb, rc, rd})
            if (v < INT_MAX) req = (req < INT_MAX) ? std::max(req, v) : v;
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return fmma(a, b, c, d, prec);
    }
    [[nodiscard]] inline Float fmms(const Float& a, const Float& b, const Float& c, const Float& d) {
        int ra = a.requestedBits(), rb = b.requestedBits();
        int rc = c.requestedBits(), rd = d.requestedBits();
        int req = INT_MAX;
        for (int v : {ra, rb, rc, rd})
            if (v < INT_MAX) req = (req < INT_MAX) ? std::max(req, v) : v;
        int prec = (req >= INT_MAX) ? Float::defaultPrecision() : Float::bitsToPrecision(req);
        return fmms(a, b, c, d, prec);
    }
    inline void sinCos(const Float& x, Float& sin_result, Float& cos_result) {
        sinCos(x, sin_result, cos_result, Float::requestedPrecision(x));
    }
    inline void sinCos(Float&& x, Float& sin_result, Float& cos_result) {
        int p = Float::requestedPrecision(x); sinCos(std::move(x), sin_result, cos_result, p);
    }
    [[nodiscard]] inline Float sinu(const Float& x, int u) { return sinu(x, u, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float sinu(Float&& x, int u)      { int p = Float::requestedPrecision(x); return sinu(std::move(x), u, p); }
    [[nodiscard]] inline Float cosu(const Float& x, int u) { return cosu(x, u, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float cosu(Float&& x, int u)      { int p = Float::requestedPrecision(x); return cosu(std::move(x), u, p); }
    [[nodiscard]] inline Float tanu(const Float& x, int u) { return tanu(x, u, Float::requestedPrecision(x)); }
    [[nodiscard]] inline Float tanu(Float&& x, int u)      { int p = Float::requestedPrecision(x); return tanu(std::move(x), u, p); }

} // namespace calx

// std 名前空間に Float の abs/sqrt + numeric_limits を追加
namespace std {
    inline calx::Float abs(const calx::Float& x) { return calx::abs(x); }
    inline calx::Float sqrt(const calx::Float& x) { return calx::sqrt(x); }

    template<>
    struct numeric_limits<calx::Float> {
        static constexpr bool is_specialized = true;
        static calx::Float epsilon() {
            // 現在のデフォルト精度に基づく machine epsilon
            int bits = calx::Float::defaultPrecision();
            // bits は十進桁数、ビットに変換: bits * log2(10) ≈ bits * 3.3219
            int prec_bits = static_cast<int>(bits * 3.3219) + 1;
            return calx::ldexp(calx::Float(1), -prec_bits);
        }
        static calx::Float min() { return calx::ldexp(calx::Float(1), -1000000); }
        static calx::Float max() { return calx::ldexp(calx::Float(1), 1000000); }
        static calx::Float lowest() { return -max(); }
        static constexpr bool is_signed = true;
        static constexpr bool is_integer = false;
        static constexpr bool is_exact = false;
        static constexpr bool has_infinity = false;
        static constexpr bool has_quiet_NaN = true;
    };
}

// --- std::formatter 特殊化 ---
#include <format>

template<>
struct std::formatter<calx::Float> {
    int precision_ = -1;
    char type_ = 'g';  // 'g'=general, 'f'=fixed, 'e'=scientific

    constexpr auto parse(std::format_parse_context& ctx) {
        auto it = ctx.begin();
        auto end = ctx.end();
        if (it != end && *it == '.') {
            ++it;
            precision_ = 0;
            while (it != end && *it >= '0' && *it <= '9') {
                precision_ = precision_ * 10 + (*it - '0');
                ++it;
            }
        }
        if (it != end && *it != '}') {
            if (*it == 'f' || *it == 'e' || *it == 'E' || *it == 'g' || *it == 'G') {
                type_ = *it;
                ++it;
            } else {
                throw std::format_error("invalid format spec for Float");
            }
        }
        return it;
    }

    auto format(const calx::Float& val, std::format_context& ctx) const {
        std::string result;
        switch (type_) {
            case 'f': result = val.toDecimalString(precision_); break;
            case 'e': case 'E': result = val.toScientificString(precision_); break;
            default: result = val.toString(precision_); break;
        }
        return std::format_to(ctx.out(), "{}", result);
    }
};

#endif // CALX_FLOAT_HPP