conv_simple for multiplication

Author: adamant-pwn

cp-algo/math/bigint.hpp

@@ -1,19 +1,19 @@
 #ifndef CP_ALGO_MATH_BIGINT_HPP
 #define CP_ALGO_MATH_BIGINT_HPP
 #include "../util/big_alloc.hpp"
-#include "../math/fft64.hpp"
+#include "../math/fft_simple.hpp"
 #include <bits/stdc++.h>
 
 namespace cp_algo::math {
     enum base_v {
-        x10 = uint64_t(1e18),
+        x10 = uint64_t(1e16),
         x16 = uint64_t(1ull << 60)
     };
     template<base_v base = x10>
     struct bigint {
-        static constexpr uint16_t digit_length = base == x10 ? 18 : 15;
+        static constexpr uint16_t digit_length = base == x10 ? 16 : 15;
         static constexpr uint16_t sub_base = base == x10 ? 10 : 16;
-        static constexpr uint32_t meta_base = base == x10 ? uint32_t(1e6) : uint32_t(1 << 20);
+        static constexpr uint32_t meta_base = base == x10 ? uint32_t(1e4) : uint32_t(1 << 15);
         big_vector<uint64_t> digits;
         bool negative;
 
@@ -122,24 +122,27 @@ namespace cp_algo::math {
         }
         void to_metabase() {
             auto N = ssize(digits);
-            digits.resize(3 * N);
+            digits.resize(4 * N);
             for (auto i = N - 1; i >= 0; i--) {
                 uint64_t val = digits[i];
-                digits[3 * i] = val % meta_base;
+                digits[4 * i] = val % meta_base;
                 val /= meta_base;
-                digits[3 * i + 1] = val % meta_base;
+                digits[4 * i + 1] = val % meta_base;
                 val /= meta_base;
-                digits[3 * i + 2] = val;
+                digits[4 * i + 2] = val % meta_base;
+                val /= meta_base;
+                digits[4 * i + 3] = val;
             }
         }
         void from_metabase() {
-            auto N = (ssize(digits) + 2) / 3;
-            digits.resize(3 * N);
+            auto N = (ssize(digits) + 3) / 4;
+            digits.resize(4 * N);
             uint64_t carry = 0;
             for (int i = 0; i < N; i++) {
-                __uint128_t val = digits[3 * i + 2];
-                val = val * meta_base + digits[3 * i + 1];
-                val = val * meta_base + digits[3 * i];
+                __uint128_t val = digits[4 * i + 3];
+                val = val * meta_base + digits[4 * i + 2];
+                val = val * meta_base + digits[4 * i + 1];
+                val = val * meta_base + digits[4 * i];
                 val += carry;
                 digits[i] = uint64_t(val % base);
                 carry = uint64_t(val / base);
@@ -172,7 +175,7 @@ namespace cp_algo::math {
             }
             to_metabase();
             other.to_metabase();
-            fft::conv64(digits, other.digits);
+            fft::conv_simple(digits, other.digits);
             from_metabase();
             return normalize();
         }

cp-algo/math/fft_simple.hpp

@@ -0,0 +1,48 @@
+#ifndef CP_ALGO_MATH_FFT_SIMPLE_HPP
+#define CP_ALGO_MATH_FFT_SIMPLE_HPP
+#include "../random/rng.hpp"
+#include "../math/common.hpp"
+#include "../math/cvector.hpp"
+CP_ALGO_SIMD_PRAGMA_PUSH
+namespace cp_algo::math::fft {
+    struct dft_simple {
+        cp_algo::math::fft::cvector cv;
+
+        dft_simple(auto const& a, size_t n): cv(n) {
+            for(size_t i = 0; i < std::min(std::size(a), n); i++) {
+                real(cv.at(i))[i % 4] = ftype(a[i]);
+                imag(cv.at(i))[i % 4] = ftype(i + n < std::size(a) ? a[i + n] : 0);
+            }
+            checkpoint("dft64 init");
+            cv.fft();
+        }
+
+        void dot(dft_simple const& t) {
+            cv.dot(t.cv);
+        }
+
+        void recover_mod(auto &res, size_t k) {
+            cv.ifft();
+            size_t n = cv.size();
+            for(size_t i = 0; i < std::min(k, n); i++) {
+                res[i] = llround(real(cv.get(i)));
+            }
+            for(size_t i = n; i < k; i++) {
+                res[i] = llround(imag(cv.get(i - n)));
+            }
+            cp_algo::checkpoint("recover mod");
+        }
+    };
+
+    // Multiplies a and b, assuming perfect precision and no overflow
+    void conv_simple(auto& a, auto const& b) {
+        size_t n = a.size(), m = b.size();
+        size_t N = std::max(flen, std::bit_ceil(n + m - 1) / 2);
+        dft_simple A(a, N), B(b, N);
+        A.dot(B);
+        a.resize(n + m - 1);
+        A.recover_mod(a, n + m - 1);
+    }
+}
+#pragma GCC pop_options
+#endif // CP_ALGO_MATH_FFT_SIMPLE_HPP

verify/bigint/hex_multiplication.test.cpp

@@ -3,7 +3,7 @@
 #pragma GCC optimize("O3,unroll-loops")
 #include <bits/allocator.h>
 #pragma GCC target("avx2")
-#define CP_ALGO_CHECKPOINT
+//#define CP_ALGO_CHECKPOINT
 #include <iostream>
 #include "blazingio/blazingio.min.hpp"
 #include "cp-algo/math/bigint.hpp"
@@ -29,5 +29,5 @@ signed main() {
     while(t--) {
         solve();
     }
-    cp_algo::checkpoint<1>();
+    //cp_algo::checkpoint<1>();
 }

verify/bigint/multiplication.test.cpp

@@ -3,7 +3,7 @@
 #pragma GCC optimize("O3,unroll-loops")
 #include <bits/allocator.h>
 #pragma GCC target("avx2")
-#define CP_ALGO_CHECKPOINT
+//#define CP_ALGO_CHECKPOINT
 #include <iostream>
 #include "blazingio/blazingio.min.hpp"
 #include "cp-algo/math/bigint.hpp"
@@ -29,5 +29,5 @@ signed main() {
     while(t--) {
         solve();
     }
-    cp_algo::checkpoint<1>();
+    //cp_algo::checkpoint<1>();
 }