GCC Code Coverage Report

Directory: ./
Coverage: low: ≥ 0% medium: ≥ 75.0% high: ≥ 90.0%
Coverage Exec / Excl / Total
Lines: 100.0% 15 / 0 / 15
Functions: 100.0% 7 / 0 / 7
Branches: -% 0 / 0 / 0
src/math/nim_prod/lib.hpp
Line Branch Exec Source
1 #pragma once
2
3 #include "../../util/alias/num/lib.hpp"
4
5 namespace tifa_libs {
6 namespace nim_prod_impl_ {
7 struct calc8 {
8 // NOLINTNEXTLINE(modernize-avoid-c-arrays)
9 u16 dp[1 << 8][1 << 8]{};
10 CEXPE calc8() NE {
11 dp[0][0] = dp[0][1] = dp[1][0] = 0, dp[1][1] = 1;
12 flt_ (u32, e, 1, 4) {
13 cu32 p = 1 << e, q = p >> 1;
14 cu16 ep = u16(1u << p), eq = u16(1u << q);
15 flt_ (u16, i, 0, ep)
16 flt_ (u16, j, i, ep) {
17 if (i < eq && j < eq) continue;
18 if (std::min(i, j) <= 1u) {
19 dp[i][j] = dp[j][i] = i * j;
20 continue;
21 }
22 cu16 iu = u16(i >> q), il = u16(i & (eq - 1)), ju = u16(j >> q), jl = u16(j & (eq - 1));
23 cu16 u = dp[iu][ju], l = dp[il][jl], ul = dp[iu ^ il][ju ^ jl], uq = dp[u][eq >> 1];
24 dp[i][j] = u16((ul ^ l) << q) ^ uq ^ l, dp[j][i] = dp[i][j];
25 }
26 }
27 }
28 } CEXP c8;
29
30 struct calc16 {
31 static CEXP u16 proot = 10279;
32 static CEXP u32 ppoly = 92191, order = 65535;
33
34 // NOLINTNEXTLINE(modernize-avoid-c-arrays)
35 u16 base[16]{}, exp[(1 << 18) + 100]{};
36 // NOLINTNEXTLINE(modernize-avoid-c-arrays)
37 u32 log[1 << 16]{};
38
39 private:
40 CEXP u16 d(u32 x) NE {
41 auto _ = x < 32768u ? 0 : ppoly;
42 return u16((x << 1) ^ _);
43 }
44 CEXP u16 naive(u16 i, u16 j) NE {
45 if (min(i, j) <= 1u) return i * j;
46 cu16 q = 8, eq = 1u << 8;
47 cu16 iu = u16(i >> q), il = u16(i & (eq - 1)), ju = u16(j >> q), jl = u16(j & (eq - 1));
48 cu16 u = c8.dp[iu][ju], l = c8.dp[il][jl], ul = c8.dp[iu ^ il][ju ^ jl], uq = c8.dp[u][eq >> 1];
49 return u16((ul ^ l) << q) ^ uq ^ l;
50 }
51
52 public:
53 CEXPE calc16() NE {
54 base[0] = 1;
55 flt_ (u32, i, 1, 16) base[i] = naive(base[i - 1], proot);
56 exp[0] = 1;
57 flt_ (u32, i, 1, order) exp[i] = d(exp[i - 1]);
58 u16* pre = exp + order + 1;
59 pre[0] = 0;
60 flt_ (u32, b, 0, 16) {
61 cu32 is = 1 << b;
62 flt_ (u32, i, is, is * 2) pre[i] = pre[i - is] ^ base[b];
63 }
64 flt_ (u32, i, 0, order) exp[i] = pre[exp[i]], log[exp[i]] = i;
65 cu32 ie_ = 2 * order + 30;
66 flt_ (u32, i, order, ie_) exp[i] = exp[i - order];
67 for (u32 i = ie_; i < sizeof(exp) / sizeof(u16); ++i) exp[i] = 0;
68 log[0] = ie_ + 1;
69 }
70 219988866 ND CEXP u16 prod(u16 i, u16 j) CNE { return exp[log[i] + log[j]]; }
71 // exp[3] = 2^{15} = 32768
72 109994433 ND CEXP u16 Hprod(u16 i, u16 j) CNE { return exp[log[i] + log[j] + 3]; }
73 36664811 ND CEXP u16 H(u16 i) CNE { return exp[log[i] + 3]; }
74 36664811 ND CEXP u16 H2(u16 i) CNE { return exp[log[i] + 6]; }
75 } CEXP c16;
76
77 CEXP u16 nimprod16(u16 i, u16 j) NE { return c16.prod(i, j); }
78 109994433 CEXP u32 nimprod32(u32 i, u32 j) NE {
79 109994433 cu16 iu = u16(i >> 16), il = u16(i & 65535), ju = u16(j >> 16), jl = u16(j & 65535);
80 109994433 cu16 l = c16.prod(il, jl), ul = c16.prod(iu ^ il, ju ^ jl), uq = c16.Hprod(iu, ju);
81 109994433 return (u32(ul ^ l) << 16) ^ uq ^ l;
82 }
83
84 // (+ : xor, x : nim product, * : integer product)
85 // i x j
86 // = (iu x ju + il x ju + iu x ji) * 2^{16}
87 // + (iu x ju x 2^{15}) + il x jl
88 // (assign ju = 2^{15}, jl = 0)
89 // = ((iu + il) x 2^{15}) * 2^{16} + (iu x 2^{15} x 2^{15})
90 36664811 CEXP u32 H(u32 i) NE {
91 36664811 cu16 iu = u16(i >> 16), il = u16(i & 65535);
92 36664811 return (u32(c16.H(iu ^ il)) << 16) ^ c16.H2(iu);
93 }
94 36664811 CEXP u64 nimprod64(u64 i, u64 j) NE {
95 36664811 cu32 iu = u32(i >> 32), il = u32(i & u32(-1)), ju = u32(j >> 32), jl = u32(j & u32(-1));
96 36664811 cu32 l = nimprod32(il, jl), ul = nimprod32(iu ^ il, ju ^ jl), uq = H(nimprod32(iu, ju));
97 36664811 return (u64(ul ^ l) << 32) ^ uq ^ l;
98 }
99 } // namespace nim_prod_impl_
100
101 using nim_prod_impl_::nimprod16, nim_prod_impl_::nimprod32, nim_prod_impl_::nimprod64;
102
103 } // namespace tifa_libs
104