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

// example_mp.cpp
// Multi-precision arithmetic demo: Int (integer), Float (floating-point), Rational (rational)
//
// Build:
//   cl /std:c++latest /EHsc /O2 /MD /utf-8 /I<calx>/include examples\example_mp.cpp
//      /link Int.lib Float.lib Rational.lib /MACHINE:X64

#include <math/core/mp/Int.hpp>
#include <math/core/mp/Int/IntGCD.hpp>
#include <math/core/mp/Int/IntPrime.hpp>
#include <math/core/mp/Int/IntCombinatorics.hpp>
#include <math/core/mp/Int/IntSqrt.hpp>
#include <math/core/mp/Float.hpp>
#include <math/core/mp/Float/FloatMath.hpp>
#include <math/core/mp/Rational.hpp>
#include <iostream>
#include <iomanip>
#include <string>
#include <chrono>

using namespace calx;

// Helper to convert Float to decimal string
// (uses integer part of pi*10^N since toDecimalString was not yet implemented)
static std::string floatToDecimal(const Float& value, int digits) {
    if (value.isNaN()) return "NaN";
    if (value.isInfinity()) return value.isNegative() ? "-Inf" : "Inf";
    if (value.isZero()) return "0.0";

    bool neg = value.isNegative();
    Float av = abs(value);

    // Compute 10^digits and multiply
    Float scale = pow(Float(10), digits, digits + 10);
    Float scaled = av * scale;
    Int int_part = scaled.toInt();
    std::string s = int_part.toString();

    // Pad with leading zeros if needed
    while (static_cast<int>(s.size()) <= digits) {
        s = "0" + s;
    }

    // Format as "X.YYYYY..."
    int int_len = static_cast<int>(s.size()) - digits;
    std::string result;
    if (neg) result += "-";
    result += s.substr(0, int_len);
    result += ".";
    result += s.substr(int_len);
    return result;
}

int main() {
    std::cout << "=== calx: Multi-Precision Arithmetic ===" << std::endl;

    // =============================================
    // Int: Multi-precision integer
    // =============================================
    std::cout << "\n--- Int: Basic Arithmetic ---" << std::endl;
    {
        Int a("123456789012345678901234567890");
        Int b("987654321098765432109876543210");
        std::cout << "a = " << a << std::endl;
        std::cout << "b = " << b << std::endl;
        std::cout << "a + b = " << a + b << std::endl;
        std::cout << "a * b = " << a * b << std::endl;
        std::cout << "b / a = " << b / a << std::endl;
        std::cout << "b % a = " << b % a << std::endl;
    }

    std::cout << "\n--- Int: Factorial ---" << std::endl;
    {
        Int n(100);
        Int result = IntCombinatorics::factorial(n);
        std::string s = result.toString();
        std::cout << "100! = " << s.substr(0, 40) << "..." << std::endl;
        std::cout << "  (" << s.size() << " digits)" << std::endl;
    }

    std::cout << "\n--- Int: GCD / LCM ---" << std::endl;
    {
        Int a(1234567890);
        Int b(9876543210LL);
        std::cout << "a = " << a << std::endl;
        std::cout << "b = " << b << std::endl;
        std::cout << "gcd(a, b) = " << IntGCD::gcd(a, b) << std::endl;
        std::cout << "lcm(a, b) = " << IntGCD::lcm(a, b) << std::endl;
    }

    std::cout << "\n--- Int: Prime Test (Miller-Rabin) ---" << std::endl;
    {
        // Mersenne prime M31 = 2^31 - 1 = 2147483647
        Int mersenne31 = pow(Int(2), Int(31)) - Int(1);
        std::cout << "2^31 - 1 = " << mersenne31 << std::endl;
        std::cout << "  isPrime: " << (IntPrime::isMillerRabinPrime(mersenne31) ? "yes" : "no")
                  << std::endl;

        // Composite number test (Int * Int)
        Int composite = Int(1234567891LL) * Int(9876543211LL);
        std::cout << composite << " isPrime: "
                  << (IntPrime::isMillerRabinPrime(composite) ? "yes" : "no") << std::endl;
    }

    std::cout << "\n--- Int: Square Root ---" << std::endl;
    {
        Int n("100000000000000000000000000000000000000");  // 10^38
        Int s = IntSqrt::sqrt(n);
        std::cout << "isqrt(10^38) = " << s << std::endl;
        std::cout << "  " << s << "^2 = " << s * s << std::endl;
    }

    // =============================================
    // Float: Multi-precision floating-point
    // =============================================
    std::cout << "\n--- Float: Mathematical Constants (50 digits) ---" << std::endl;
    {
        int prec = 60;  // compute at 60-digit precision with margin
        int show = 50;  // display 50 digits
        Float pi_val = Float::pi(prec);
        Float e_val = Float::e(prec);
        Float log2_val = Float::log2(prec);
        Float euler_val = Float::euler(prec);
        Float sqrt2_val = Float::sqrt2(prec);

        std::cout << "pi    = " << floatToDecimal(pi_val, show) << std::endl;
        std::cout << "e     = " << floatToDecimal(e_val, show) << std::endl;
        std::cout << "log2  = " << floatToDecimal(log2_val, show) << std::endl;
        std::cout << "gamma = " << floatToDecimal(euler_val, show) << std::endl;
        std::cout << "sqrt2 = " << floatToDecimal(sqrt2_val, show) << std::endl;
    }

    std::cout << "\n--- Float: Pi Calculation (Chudnovsky, 10000 digits) ---" << std::endl;
    {
        int prec = 10010;
        auto start = std::chrono::high_resolution_clock::now();
        Float pi_val = Float::pi(prec);
        auto end = std::chrono::high_resolution_clock::now();
        auto ms = std::chrono::duration_cast<std::chrono::milliseconds>(end - start).count();

        std::string pi_str = floatToDecimal(pi_val, 10000);
        std::cout << "pi (first 200 of " << pi_str.size() - 2 << " digits):" << std::endl;
        for (size_t i = 0; i < 200 && i < pi_str.size(); i += 80) {
            size_t len = std::min<size_t>(80, pi_str.size() - i);
            if (i + len > 200) len = 200 - i;
            std::cout << "  " << pi_str.substr(i, len) << std::endl;
        }
        std::cout << "  ... (" << pi_str.size() - 2 << " digits total)" << std::endl;
        std::cout << "  Computation time: " << ms << " ms" << std::endl;
    }

    std::cout << "\n--- Float: Math Functions ---" << std::endl;
    {
        int prec = 50;
        std::cout << std::setprecision(15);
        std::cout << "exp(1)    = " << exp(Float(1), prec).toDouble()
                  << "  (exact: 2.71828182845905)" << std::endl;
        std::cout << "log(2)    = " << log(Float(2), prec).toDouble()
                  << "  (exact: 0.693147180559945)" << std::endl;
        std::cout << "sin(pi/6) = " << sin(Float::pi(prec) / Float(6), prec).toDouble()
                  << "  (exact: 0.5)" << std::endl;
        std::cout << "sqrt(2)   = " << sqrt(Float(2), prec).toDouble()
                  << "  (exact: 1.41421356237310)" << std::endl;
    }

    // =============================================
    // Rational: Exact rational arithmetic
    // =============================================
    std::cout << "\n--- Rational: Exact Arithmetic ---" << std::endl;
    {
        Rational a(1, 3);   // 1/3
        Rational b(1, 6);   // 1/6
        std::cout << "a = " << a << std::endl;
        std::cout << "b = " << b << std::endl;
        std::cout << "a + b = " << a + b << std::endl;   // 1/2
        std::cout << "a * b = " << a * b << std::endl;   // 1/18
        std::cout << "a / b = " << a / b << std::endl;   // 2
    }

    std::cout << "\n--- Rational: Harmonic Series ---" << std::endl;
    {
        // Exact computation of harmonic series H_20 = 1 + 1/2 + 1/3 + ... + 1/20
        Rational sum;
        for (int i = 1; i <= 20; ++i) {
            sum += Rational(1, i);
        }
        std::cout << "H_20 = " << sum << std::endl;
        std::cout << "     = " << sum.toDecimal(30) << "..." << std::endl;
    }

    std::cout << "\n--- Rational: Continued Fraction ---" << std::endl;
    {
        // 355/113 (an excellent rational approximation of pi)
        Rational pi_approx(355, 113);
        auto cf = pi_approx.toContinuedFraction();
        std::cout << "355/113 = " << pi_approx.toDecimal(20) << std::endl;
        std::cout << "  Continued fraction: [";
        for (size_t i = 0; i < cf.size(); ++i) {
            if (i > 0) std::cout << "; ";
            std::cout << cf[i];
        }
        std::cout << "]" << std::endl;
    }

    return 0;
}