mode math algorithms to math folder

math/fibonacci_fast.cpp

@@ -0,0 +1,52 @@
+/**
+ * @file
+ * @brief Faster computation of Fibonacci series
+ *
+ * An efficient way to calculate nth fibonacci number faster and simpler than
+ * \f$O(n\log n)\f$ method of matrix exponentiation This works by using both
+ * recursion and dynamic programming. as 93rd fibonacci exceeds 19 digits, which
+ * cannot be stored in a single long long variable, we can only use it till 92nd
+ * fibonacci we can use it for 10000th fibonacci etc, if we implement
+ * bigintegers. This algorithm works with the fact that nth fibonacci can easily
+ * found if we have already found n/2th or (n+1)/2th fibonacci It is a property
+ * of fibonacci similar to matrix exponentiation.
+ *
+ * @see fibonacci_large.cpp, fibonacci.cpp
+ */
+
+#include <cinttypes>
+#include <cstdio>
+#include <iostream>
+
+/** maximum number that can be computed - The result after 93 cannot be stored
+ * in a `uint64_t` data type. */
+const uint64_t MAX = 93;
+
+/** Array of computed fibonacci numbers */
+uint64_t f[MAX] = {0};
+
+/** Algorithm */
+uint64_t fib(uint64_t n) {
+    if (n == 0)
+        return 0;
+    if (n == 1 || n == 2)
+        return (f[n] = 1);
+
+    if (f[n])
+        return f[n];
+
+    uint64_t k = (n % 2 != 0) ? (n + 1) / 2 : n / 2;
+
+    f[n] = (n % 2 != 0) ? (fib(k) * fib(k) + fib(k - 1) * fib(k - 1))
+                        : (2 * fib(k - 1) + fib(k)) * fib(k);
+    return f[n];
+}
+
+/** Main function */
+int main() {
+    // Main Function
+    for (uint64_t i = 1; i < 93; i++) {
+        std::cout << i << " th fibonacci number is " << fib(i) << std::endl;
+    }
+    return 0;
+}

math/fibonacci_large.cpp

@@ -0,0 +1,84 @@
+/**
+ * @file
+ * @brief Computes N^th Fibonacci number given as
+ * input argument. Uses custom build arbitrary integers library
+ * to perform additions and other operations.
+ *
+ * Took 0.608246 seconds to compute 50,000^th Fibonacci
+ * number that contains 10450 digits!
+ *
+ * @see fibonacci.cpp, fibonacci_fast.cpp
+ */
+
+#include <cinttypes>
+#include <ctime>
+#include <iostream>
+
+#include "./large_number.h"
+
+/** Compute fibonacci numbers using the relation
+ * \f[f(n)=f(n-1)+f(n-2)\f]
+ * and returns the result as a large_number type.
+ */
+large_number fib(uint64_t n) {
+    large_number f0(1);
+    large_number f1(1);
+
+    do {
+        large_number f2 = f1;
+        f1 += f0;
+        f0 = f2;
+        n--;
+    } while (n > 2);  // since we start from 2
+
+    return f1;
+}
+
+int main(int argc, char *argv[]) {
+    uint64_t N;
+    if (argc == 2) {
+        N = strtoull(argv[1], NULL, 10);
+    } else {
+        std::cout << "Enter N: ";
+        std::cin >> N;
+    }
+
+    clock_t start_time = std::clock();
+    large_number result = fib(N);
+    clock_t end_time = std::clock();
+    double time_taken = static_cast<double>(end_time - start_time) /
+                        static_cast<double>(CLOCKS_PER_SEC);
+
+    std::cout << std::endl
+              << N << "^th Fibonacci number: " << result << std::endl
+              << "Number of digits: " << result.num_digits() << std::endl
+              << "Time taken: " << std::scientific << time_taken << " s"
+              << std::endl;
+
+    N = 5000;
+    if (fib(N) ==
+        large_number(
+            "387896845438832563370191630832590531208212771464624510616059721489"
+            "555013904403709701082291646221066947929345285888297381348310200895"
+            "498294036143015691147893836421656394410691021450563413370655865623"
+            "825465670071252592990385493381392883637834751890876297071203333705"
+            "292310769300851809384980180384781399674888176555465378829164426891"
+            "298038461377896902150229308247566634622492307188332480328037503913"
+            "035290330450584270114763524227021093463769910400671417488329842289"
+            "149127310405432875329804427367682297724498774987455569190770388063"
+            "704683279481135897373999311010621930814901857081539785437919530561"
+            "751076105307568878376603366735544525884488624161921055345749367589"
+            "784902798823435102359984466393485325641195222185956306047536464547"
+            "076033090242080638258492915645287629157575914234380914230291749108"
+            "898415520985443248659407979357131684169286803954530954538869811466"
+            "508206686289742063932343848846524098874239587380197699382031717420"
+            "893226546887936400263079778005875912967138963421425257911687275560"
+            "0360311370547754724604639987588046985178408674382863125"))
+        std::cout << "Test for " << N << "^th Fibonacci number passed!"
+                  << std::endl;
+    else
+        std::cerr << "Test for " << N << "^th Fibonacci number failed!"
+                  << std::endl;
+
+    return 0;
+}

math/large_factorial.cpp

@@ -0,0 +1,117 @@
+/**
+ * @file
+ * @brief Compute factorial of any arbitratily large number/
+ *
+ * @see factorial.cpp
+ */
+#include <cstring>
+#include <ctime>
+#include <iostream>
+
+#include "./large_number.h"
+
+/** Test implementation for 10! Result must be 3628800.
+ * @returns True if test pass else False
+ */
+bool test1() {
+    std::cout << "---- Check 1\t";
+    unsigned int i, number = 10;
+    large_number result;
+    for (i = 2; i <= number; i++) /* Multiply every number from 2 thru N */
+        result *= i;
+
+    const char *known_reslt = "3628800";
+
+    /* check 1 */
+    if (strlen(known_reslt) != result.num_digits()) {
+        std::cerr << "Result lengths dont match! " << strlen(known_reslt)
+                  << " != " << result.num_digits() << std::endl;
+        return false;
+    }
+
+    const size_t N = result.num_digits();
+    for (i = 0; i < N; i++) {
+        if (known_reslt[i] != result.digit_char(i)) {
+            std::cerr << i << "^th digit mismatch! " << known_reslt[i]
+                      << " != " << result.digit_char(i) << std::endl;
+            return false;
+        }
+    }
+
+    std::cout << "Passed!" << std::endl;
+    return true;
+}
+
+/** Test implementation for 100! The result is the 156 digit number:
+ * ```
+ * 9332621544394415268169923885626670049071596826438162146859296389521759
+ * 9993229915608941463976156518286253697920827223758251185210916864000000
+ * 000000000000000000
+ * ```
+ * @returns True if test pass else False
+ */
+bool test2() {
+    std::cout << "---- Check 2\t";
+    unsigned int i, number = 100;
+    large_number result;
+    for (i = 2; i <= number; i++) /* Multiply every number from 2 thru N */
+        result *= i;
+
+    const char *known_reslt =
+        "9332621544394415268169923885626670049071596826438162146859296389521759"
+        "9993229915608941463976156518286253697920827223758251185210916864000000"
+        "000000000000000000";
+
+    /* check 1 */
+    if (strlen(known_reslt) != result.num_digits()) {
+        std::cerr << "Result lengths dont match! " << strlen(known_reslt)
+                  << " != " << result.num_digits() << std::endl;
+        return false;
+    }
+
+    const size_t N = result.num_digits();
+    for (i = 0; i < N; i++) {
+        if (known_reslt[i] != result.digit_char(i)) {
+            std::cerr << i << "^th digit mismatch! " << known_reslt[i]
+                      << " != " << result.digit_char(i) << std::endl;
+            return false;
+        }
+    }
+
+    std::cout << "Passed!" << std::endl;
+    return true;
+}
+
+/**
+ * Main program
+ **/
+int main(int argc, char *argv[]) {
+    int number, i;
+
+    if (argc == 2) {
+        number = atoi(argv[1]);
+    } else {
+        std::cout << "Enter the value of n(n starts from 0 ): ";
+        std::cin >> number;
+    }
+
+    large_number result;
+
+    std::clock_t start_time = std::clock();
+    for (i = 2; i <= number; i++) /* Multiply every number from 2 thru N */
+        result *= i;
+    std::clock_t end_time = std::clock();
+    double time_taken =
+        static_cast<double>(end_time - start_time) / CLOCKS_PER_SEC;
+
+    std::cout << number << "! = " << result << std::endl
+              << "Number of digits: " << result.num_digits() << std::endl
+              << "Time taken: " << std::scientific << time_taken << " s"
+              << std::endl;
+
+    test1();
+    test2();
+    result.test();
+
+    return 0;
+}

math/large_number.h

@@ -0,0 +1,287 @@
+/**
+ * @file
+ * @brief Library to perform arithmatic operations on arbitrarily large
+ * numbers.
+ */
+
+#ifndef OTHERS_LARGE_NUMBER_H_
+#define OTHERS_LARGE_NUMBER_H_
+#include <algorithm>
+#include <cassert>
+#include <cinttypes>
+#include <cstring>
+#include <iostream>
+#include <type_traits>
+#include <vector>
+
+/**
+ * Store large unsigned numbers as a C++ vector
+ * The class provides convenience functions to add a
+ * digit to the number, perform multiplication of
+ * large number with long unsigned integers.
+ **/
+class large_number {
+ public:
+    /**< initializer with value = 1 */
+    large_number() { _digits.push_back(1); }
+
+    // /**< initializer from an integer */
+    // explicit large_number(uint64_t n) {
+    //     uint64_t carry = n;
+    //     do {
+    //         add_digit(carry % 10);
+    //         carry /= 10;
+    //     } while (carry != 0);
+    // }
+
+    /**< initializer from an integer */
+    explicit large_number(int n) {
+        int carry = n;
+        do {
+            add_digit(carry % 10);
+            carry /= 10;
+        } while (carry != 0);
+    }
+
+    /**< initializer from another large_number */
+    large_number(const large_number &a) : _digits(a._digits) {}
+
+    /**< initializer from a vector */
+    explicit large_number(std::vector<unsigned char> &vec) : _digits(vec) {}
+
+    /**< initializer from a string */
+    explicit large_number(char const *number_str) {
+        for (size_t i = strlen(number_str); i > 0; i--) {
+            unsigned char a = number_str[i - 1] - '0';
+            if (a >= 0 && a <= 9)
+                _digits.push_back(a);
+        }
+    }
+
+    /**
+     * Function to check implementation
+     **/
+    static bool test() {
+        std::cout << "------ Checking `large_number` class implementations\t"
+                  << std::endl;
+        large_number a(40);
+        // 1. test multiplication
+        a *= 10;
+        if (a != large_number(400)) {
+            std::cerr << "\tFailed 1/6 (" << a << "!=400)" << std::endl;
+            return false;
+        }
+        std::cout << "\tPassed 1/6...";
+        // 2. test compound addition with integer
+        a += 120;
+        if (a != large_number(520)) {
+            std::cerr << "\tFailed 2/6 (" << a << "!=520)" << std::endl;
+            return false;
+        }
+        std::cout << "\tPassed 2/6...";
+        // 3. test compound multiplication again
+        a *= 10;
+        if (a != large_number(5200)) {
+            std::cerr << "\tFailed 3/6 (" << a << "!=5200)" << std::endl;
+            return false;
+        }
+        std::cout << "\tPassed 3/6...";
+        // 4. test increment (prefix)
+        ++a;
+        if (a != large_number(5201)) {
+            std::cerr << "\tFailed 4/6 (" << a << "!=5201)" << std::endl;
+            return false;
+        }
+        std::cout << "\tPassed 4/6...";
+        // 5. test increment (postfix)
+        a++;
+        if (a != large_number(5202)) {
+            std::cerr << "\tFailed 5/6 (" << a << "!=5202)" << std::endl;
+            return false;
+        }
+        std::cout << "\tPassed 5/6...";
+        // 6. test addition with another large number
+        a = a + large_number("7000000000000000000000000000000");
+        if (a != large_number("7000000000000000000000000005202")) {
+            std::cerr << "\tFailed 6/6 (" << a
+                      << "!=7000000000000000000000000005202)" << std::endl;
+            return false;
+        }
+        std::cout << "\tPassed 6/6..." << std::endl;
+        return true;
+    }
+
+    /**
+     * add a digit at MSB to the large number
+     **/
+    void add_digit(unsigned int value) {
+        if (value > 9) {
+            std::cerr << "digit > 9!!\n";
+            exit(EXIT_FAILURE);
+        }
+
+        _digits.push_back(value);
+    }
+
+    /**
+     * Get number of digits in the number
+     **/
+    const size_t num_digits() const { return _digits.size(); }
+
+    /**
+     * operator over load to access the
+     * i^th digit conveniently and also
+     * assign value to it
+     **/
+    inline unsigned char &operator[](size_t n) { return this->_digits[n]; }
+
+    inline const unsigned char &operator[](size_t n) const {
+        return this->_digits[n];
+    }
+
+    /**
+     * operator overload to compare two numbers
+     **/
+    friend std::ostream &operator<<(std::ostream &out, const large_number &a) {
+        for (size_t i = a.num_digits(); i > 0; i--)
+            out << static_cast<int>(a[i - 1]);
+        return out;
+    }
+
+    /**
+     * operator overload to compare two numbers
+     **/
+    friend bool operator==(large_number const &a, large_number const &b) {
+        size_t N = a.num_digits();
+        if (N != b.num_digits())
+            return false;
+        for (size_t i = 0; i < N; i++)
+            if (a[i] != b[i])
+                return false;
+        return true;
+    }
+
+    /**
+     * operator overload to compare two numbers
+     **/
+    friend bool operator!=(large_number const &a, large_number const &b) {
+        return !(a == b);
+    }
+
+    /**
+     * operator overload to increment (prefix)
+     **/
+    large_number &operator++() {
+        (*this) += 1;
+        return *this;
+    }
+
+    /**
+     * operator overload to increment (postfix)
+     **/
+    large_number &operator++(int) {
+        static large_number tmp(_digits);
+        ++(*this);
+        return tmp;
+    }
+
+    /**
+     * operator overload to add
+     **/
+    large_number &operator+=(large_number n) {
+        // if adding with another large_number
+        large_number *b = reinterpret_cast<large_number *>(&n);
+        const size_t max_L = std::max(this->num_digits(), b->num_digits());
+        unsigned int carry = 0;
+        size_t i;
+        for (i = 0; i < max_L || carry != 0; i++) {
+            if (i < b->num_digits())
+                carry += (*b)[i];
+            if (i < this->num_digits())
+                carry += (*this)[i];
+            if (i < this->num_digits())
+                (*this)[i] = carry % 10;
+            else
+                this->add_digit(carry % 10);
+            carry /= 10;
+        }
+        return *this;
+    }
+
+    large_number &operator+=(int n) { return (*this) += large_number(n); }
+    // large_number &operator+=(uint64_t n) { return (*this) += large_number(n);
+    // }
+
+    /**
+     * operator overload to perform addition
+     **/
+    template <class T>
+    friend large_number &operator+(const large_number &a, const T &b) {
+        static large_number c = a;
+        c += b;
+        return c;
+    }
+
+    /**
+     * assignment operator
+     **/
+    large_number &operator=(const large_number &b) {
+        this->_digits = b._digits;
+        return *this;
+    }
+
+    /**
+     * operator overload to increment
+     **/
+    template <class T>
+    large_number &operator*=(const T n) {
+        static_assert(std::is_integral<T>::value,
+                      "Must be integer addition unsigned integer types.");
+        this->multiply(n);
+        return *this;
+    }
+
+    /**
+     * returns i^th digit as an ASCII character
+     **/
+    const char digit_char(size_t i) const {
+        return _digits[num_digits() - i - 1] + '0';
+    }
+
+ private:
+    /**
+     * multiply large number with another integer and
+     * store the result in the same large number
+     **/
+    template <class T>
+    void multiply(const T n) {
+        static_assert(std::is_integral<T>::value,
+                      "Can only have integer types.");
+        // assert(!(std::is_signed<T>::value)); //, "Implemented only for
+        // unsigned integer types.");
+
+        size_t i;
+        uint64_t carry = 0, temp;
+        for (i = 0; i < this->num_digits(); i++) {
+            temp = (*this)[i] * n;
+            temp += carry;
+            if (temp < 10) {
+                carry = 0;
+            } else {
+                carry = temp / 10;
+                temp = temp % 10;
+            }
+            (*this)[i] = temp;
+        }
+
+        while (carry != 0) {
+            this->add_digit(carry % 10);
+            carry /= 10;
+        }
+    }
+
+    std::vector<unsigned char>
+        _digits; /**< where individual digits are stored */
+};
+
+#endif  // OTHERS_LARGE_NUMBER_H_