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

// UInt128.hpp
// Definition of the 128-bit unsigned integer type
//
// This file provides the 128-bit unsigned integer type that was split out from basic_types.hpp.

#ifndef SANGI_UINT128_HPP
#define SANGI_UINT128_HPP

#include <cstdint>
#include <utility>
#include <stdexcept>
#include <bit>
#include <math/core/mp/Int/compiler_intrinsics.hpp>

namespace sangi {

    /**
     * @brief 128-bit unsigned integer type (compiler independent)
     *
     * This struct provides a 128-bit unsigned integer type usable on any compiler.
     * It can be used as an alternative to GCC/Clang's __uint128_t.
     */
    struct UInt128 {
        uint64_t low;   ///< Lower 64 bits
        uint64_t high;  ///< Upper 64 bits

        UInt128() : low(0), high(0) {}
        UInt128(uint64_t l) : low(l), high(0) {}
        UInt128(uint64_t h, uint64_t l) : low(l), high(h) {}

        /**
         * @brief Construct a UInt128 from upper 64 bits
         * @param high Upper 64-bit value
         * @return The constructed UInt128 value
         */
        static UInt128 shift64(uint64_t high) {
            return UInt128(high, 0);
        }

        /**
         * @brief Addition
         * @param other Value to add
         * @return Reference to the addition result
         */
        UInt128& operator+=(const UInt128& other) {
            uint64_t old_low = low;
            low += other.low;
            // Carry calculation: a carry occurs if the post-addition value is smaller than the pre-addition value
            high += other.high + (low < old_low ? 1 : 0);
            return *this;
        }

        /**
         * @brief Addition with a uint64_t value
         * @param value Value to add
         * @return Reference to the addition result
         */
        UInt128& operator+=(uint64_t value) {
            uint64_t old_low = low;
            low += value;
            // Carry calculation: a carry occurs if the post-addition value is smaller than the pre-addition value
            if (low < old_low) {
                high++;
            }
            return *this;
        }

        /**
         * @brief Subtraction
         * @param other Value to subtract
         * @return Reference to the subtraction result
         */
        UInt128& operator-=(const UInt128& other) {
            uint64_t old_low = low;
            low -= other.low;
            // Borrow calculation: a borrow occurs if the post-subtraction value is larger than the pre-subtraction value
            high -= other.high + (low > old_low ? 1 : 0);
            return *this;
        }

        /**
         * @brief Subtraction with a uint64_t value
         * @param value Value to subtract
         * @return Reference to the subtraction result
         */
        UInt128& operator-=(uint64_t value) {
            uint64_t old_low = low;
            low -= value;
            // Borrow calculation: a borrow occurs if the post-subtraction value is larger than the pre-subtraction value
            if (low > old_low) {
                high--;
            }
            return *this;
        }

        /**
         * @brief Multiplication of two 64-bit integers
         * @param a First operand
         * @param b Second operand
         * @return 128-bit result of a * b
         */
        /**
         * @brief 64-bit x 64-bit -> 128-bit multiplication
         *
         * MSVC: use _umul128 intrinsic (one instruction)
         * GCC/Clang: use __uint128_t or __builtin_umul_overflow
         * Fallback: split into 32-bit words (legacy implementation)
         */
        static UInt128 multiply(uint64_t a, uint64_t b) {
#ifdef _MSC_VER
            // MSVC: use _umul128 intrinsic
            uint64_t high;
            uint64_t low = _umul128(a, b, &high);
            return UInt128(high, low);
#elif defined(__SIZEOF_INT128__)
            // GCC/Clang: use __uint128_t
            __uint128_t result = static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
            return UInt128(static_cast<uint64_t>(result >> 64), static_cast<uint64_t>(result));
#else
            // Fallback: software implementation (split into 32-bit words)
            constexpr uint64_t MASK32 = 0xFFFFFFFF;

            uint64_t a_low = a & MASK32;
            uint64_t a_high = a >> 32;
            uint64_t b_low = b & MASK32;
            uint64_t b_high = b >> 32;

            uint64_t ll = a_low * b_low;
            uint64_t lh = a_low * b_high;
            uint64_t hl = a_high * b_low;
            uint64_t hh = a_high * b_high;

            // mid = ll>>32 + (lh & MASK32) + (hl & MASK32)
            uint64_t mid = (ll >> 32) + (lh & MASK32) + (hl & MASK32);
            uint64_t low = (ll & MASK32) | (mid << 32);

            // Build high including the carry
            uint64_t high = hh + (lh >> 32) + (hl >> 32) + (mid >> 32);

            return UInt128(high, low);
#endif
        }

        /**
         * @brief Division (128 bits / 64 bits)
         * @param dividend Dividend
         * @param divisor Divisor
         * @return Pair of quotient and remainder
         */
        static std::pair<uint64_t, uint64_t> divide(const UInt128& dividend, uint64_t divisor) {
            // Throw an exception when divisor is 0
            if (divisor == 0) {
                throw std::invalid_argument("Division by zero");
            }

            // Overflow when the high word of dividend is >= divisor
            if (dividend.high >= divisor) {
                throw std::overflow_error("Division overflow");
            }

            // Standard division processing
            uint64_t q = 0;
            uint64_t r = 0;

            // Process from the high digits
            for (int i = 63; i >= 0; --i) {
                // Left-shift the remainder
                r = (r << 1) | ((dividend.high >> i) & 1);

                // Compute the quotient
                if (r >= divisor) {
                    r -= divisor;
                    q |= (1ULL << i);
                }
            }

            // Process the lower digits
            for (int i = 63; i >= 0; --i) {
                // Left-shift the remainder
                r = (r << 1) | ((dividend.low >> i) & 1);

                // Compute the quotient
                if (r >= divisor) {
                    r -= divisor;
                    q |= (1ULL << i);
                }
            }

            return std::make_pair(q, r); // Quotient and remainder
        }

        /**
         * @brief Fast algorithm: 128 bits / 64 bits = 64-bit quotient and 64-bit remainder
         *
         * MSVC: use _udiv128 intrinsic (one instruction)
         * GCC/Clang: __uint128_t or manual implementation
         * Fallback: bit-by-bit trial subtraction (legacy implementation)
         *
         * @param high Upper 64 bits of the dividend
         * @param low Lower 64 bits of the dividend
         * @param divisor Divisor
         * @return Pair of quotient and remainder
         */
        static std::pair<uint64_t, uint64_t> divmod_fast(uint64_t high, uint64_t low, uint64_t divisor) {
            // Throw an exception when divisor is 0
            if (divisor == 0) {
                throw std::invalid_argument("Division by zero");
            }

            // Overflow when high >= divisor
            if (high >= divisor) {
                throw std::overflow_error("Division overflow");
            }

#if defined(_MSC_VER) && defined(_M_X64)
            // MSVC x64: use _udiv128 intrinsic
            uint64_t remainder;
            uint64_t quotient = _udiv128(high, low, divisor, &remainder);
            return std::make_pair(quotient, remainder);
#elif defined(__SIZEOF_INT128__)
            // GCC/Clang: use __uint128_t
            __uint128_t dividend = (static_cast<__uint128_t>(high) << 64) | low;
            uint64_t quotient = static_cast<uint64_t>(dividend / divisor);
            uint64_t remainder = static_cast<uint64_t>(dividend % divisor);
            return std::make_pair(quotient, remainder);
#else
            // Fallback: bit-by-bit trial subtraction (software implementation)
            uint64_t quotient = 0;
            uint64_t remainder = high;

            // Process the upper 64 bits
            for (int i = 0; i < 64; ++i) {
                // When the MSB of remainder is set, left-shifting overflows.
                // In that case, the true post-shift value is >= 2^64, hence always >= divisor.
                bool overflow = (remainder >> 63) & 1;

                // Left-shift the remainder and bring in the high bit of low
                remainder = (remainder << 1) | ((low >> 63) & 1);
                low <<= 1;

                // Compute the quotient bit
                quotient <<= 1;
                if (overflow || remainder >= divisor) {
                    remainder -= divisor;
                    quotient |= 1;
                }
            }

            return std::make_pair(quotient, remainder);
#endif
        }

        /**
         * @brief Comparison (equal)
         * @param other Comparison target
         * @return true if equal
         */
        bool operator==(const UInt128& other) const {
            return low == other.low && high == other.high;
        }

        /**
         * @brief Comparison (not equal)
         * @param other Comparison target
         * @return true if not equal
         */
        bool operator!=(const UInt128& other) const {
            return !(*this == other);
        }

        /**
         * @brief Comparison (ordering)
         * @param other Comparison target
         * @return true if greater
         */
        bool operator>(const UInt128& other) const {
            return high > other.high || (high == other.high && low > other.low);
        }

        /**
         * @brief Comparison (ordering)
         * @param other Comparison target
         * @return true if less
         */
        bool operator<(const UInt128& other) const {
            return high < other.high || (high == other.high && low < other.low);
        }

        /**
         * @brief Comparison (>=)
         * @param other Comparison target
         * @return true if greater than or equal
         */
        bool operator>=(const UInt128& other) const {
            return !(*this < other);
        }

        /**
         * @brief Comparison (<=)
         * @param other Comparison target
         * @return true if less than or equal
         */
        bool operator<=(const UInt128& other) const {
            return !(*this > other);
        }

        /**
         * @brief Comparison with uint64_t (equal)
         * @param value Comparison target
         * @return true if equal
         */
        bool operator==(uint64_t value) const {
            return high == 0 && low == value;
        }

        /**
         * @brief Comparison with uint64_t (not equal)
         * @param value Comparison target
         * @return true if not equal
         */
        bool operator!=(uint64_t value) const {
            return !(*this == value);
        }

        /**
         * @brief Comparison with uint64_t (ordering)
         * @param value Comparison target
         * @return true if greater
         */
        bool operator>(uint64_t value) const {
            return high > 0 || low > value;
        }

        /**
         * @brief Comparison with uint64_t (ordering)
         * @param value Comparison target
         * @return true if less
         */
        bool operator<(uint64_t value) const {
            return high == 0 && low < value;
        }

        /**
         * @brief Comparison with uint64_t (>=)
         * @param value Comparison target
         * @return true if greater than or equal
         */
        bool operator>=(uint64_t value) const {
            return !(*this < value);
        }

        /**
         * @brief Comparison with uint64_t (<=)
         * @param value Comparison target
         * @return true if less than or equal
         */
        bool operator<=(uint64_t value) const {
            return !(*this > value);
        }
    };

    // Non-member operator overloads

    /**
     * @brief Addition
     * @param lhs Left-hand side
     * @param rhs Right-hand side
     * @return Addition result
     */
    inline UInt128 operator+(const UInt128& lhs, const UInt128& rhs) {
        UInt128 result = lhs;
        result += rhs;
        return result;
    }

    /**
     * @brief Addition of UInt128 and uint64_t
     * @param lhs Left-hand side
     * @param rhs Right-hand side
     * @return Addition result
     */
    inline UInt128 operator+(const UInt128& lhs, uint64_t rhs) {
        UInt128 result = lhs;
        result += rhs;
        return result;
    }

    /**
     * @brief Addition of uint64_t and UInt128
     * @param lhs Left-hand side
     * @param rhs Right-hand side
     * @return Addition result
     */
    inline UInt128 operator+(uint64_t lhs, const UInt128& rhs) {
        UInt128 result(0, lhs);
        result += rhs;
        return result;
    }

    /**
     * @brief Subtraction
     * @param lhs Left-hand side
     * @param rhs Right-hand side
     * @return Subtraction result
     */
    inline UInt128 operator-(const UInt128& lhs, const UInt128& rhs) {
        UInt128 result = lhs;
        result -= rhs;
        return result;
    }

    /**
     * @brief Subtraction of UInt128 and uint64_t
     * @param lhs Left-hand side
     * @param rhs Right-hand side
     * @return Subtraction result
     */
    inline UInt128 operator-(const UInt128& lhs, uint64_t rhs) {
        UInt128 result = lhs;
        result -= rhs;
        return result;
    }

    /**
     * @brief Subtraction of uint64_t and UInt128
     * @param lhs Left-hand side
     * @param rhs Right-hand side
     * @return Subtraction result
     */
    inline UInt128 operator-(uint64_t lhs, const UInt128& rhs) {
        if (rhs.high > 0 || lhs < rhs.low) {
            throw std::underflow_error("Unsigned subtraction underflow");
        }
        return UInt128(0, lhs - rhs.low);
    }

} // namespace sangi

#endif // SANGI_UINT128_HPP