GCC Code Coverage Report

Directory: src/
Coverage: low: ≥ 0% medium: ≥ 75.0% high: ≥ 90.0%
Coverage Exec / Excl / Total
Lines: 76.7% 178 / 0 / 232
Functions: 90.4% 66 / 0 / 73
Branches: 86.0% 294 / 35 / 377
num/modnum.hpp
Line Branch Exec Source
1 #pragma once
2
3 #include <cassert>
4 #include <istream>
5 #include <ostream>
6 #include <cstdint>
7 #include <concepts>
8 #include <limits>
9 #include <type_traits>
10
11 namespace wala {
12
13 10188753 template <typename T> T mod_inv_in_range(T a, T m) {
14 // assert(0 <= a && a < m);
15 10188753 T x = a, y = m;
16 // abs coeff of a in x and y (they're always opposite sign)
17 10188753 T vx = 1, vy = 0;
18 10188753 bool swap = false;
19
8/8
unsigned int wala::mod_inv_in_range<unsigned int>(unsigned int, unsigned int):
✓ Branch 4 → 3 taken 30947925 times.
✓ Branch 4 → 5 taken 10162902 times.
✓ Branch 10 → 7 taken 82025 times.
✓ Branch 10 → 11 taken 16125 times.
unsigned long wala::mod_inv_in_range<unsigned long>(unsigned long, unsigned long):
✓ Branch 4 → 3 taken 14270 times.
✓ Branch 4 → 5 taken 628 times.
✓ Branch 10 → 7 taken 64324 times.
✓ Branch 10 → 11 taken 9098 times.
 41297297 while (x) {
20 31108544 T k = y / x;
21 31108544 y %= x;
22 31108544 vy += k * vx;
23 31108544 std::swap(x, y);
24 31108544 std::swap(vx, vy);
25 31108544 swap ^= 1;
26 }
27 10188753 assert(y == 1);
28
8/8
unsigned int wala::mod_inv_in_range<unsigned int>(unsigned int, unsigned int):
✓ Branch 7 → 8 taken 9495547 times.
✓ Branch 7 → 9 taken 667355 times.
✓ Branch 13 → 14 taken 5501 times.
✓ Branch 13 → 15 taken 10624 times.
unsigned long wala::mod_inv_in_range<unsigned long>(unsigned long, unsigned long):
✓ Branch 7 → 8 taken 476 times.
✓ Branch 7 → 9 taken 152 times.
✓ Branch 13 → 14 taken 5432 times.
✓ Branch 13 → 15 taken 3666 times.
 10188753 return swap ? vy : m - vy;
29 }
30
31 template <typename T> struct extended_gcd_result {
32 T gcd;
33 T coeff_a, coeff_b;
34 };
35 2018 template <typename T> extended_gcd_result<T> extended_gcd(T a, T b) {
36 2018 T x = a, y = b;
37 // coeff of a and b in x and y
38 2018 T ax = 1, ay = 0;
39 2018 T bx = 0, by = 1;
40
2/2
✓ Branch 14 → 9 taken 6604 times.
✓ Branch 14 → 15 taken 2018 times.
 8622 while (x) {
41
1/2
✗ Branch 9 → 10 not taken.
✓ Branch 9 → 11 taken 6604 times.
 6604 T k = y / x;
42 6604 y %= x;
43 6604 ay -= k * ax;
44 6604 by -= k * bx;
45 6604 std::swap(x, y);
46 6604 std::swap(ax, ay);
47 6604 std::swap(bx, by);
48 }
49 2018 return {y, ay, by};
50 }
51
52 ✗ template <typename T> T mod_inv(T a, T m) {
53 ✗ a %= m;
54 ✗ a = a < 0 ? a + m : a;
55 ✗ return mod_inv_in_range(a, m);
56 }
57
58 // Derives the boilerplate operator surface of a number type from its compound
59 // ops, ==, neg(), and inv().
60 // Bodies are only instantiated on use, so a type may omit some of the
61 // underlying pieces if the corresponding derived ops are never called.
62 template <typename Self>
63 struct num_ops {
64 ✗ Self operator+ () const { return static_cast<const Self&>(*this); }
65 32392196 Self operator- () const { return static_cast<const Self&>(*this).neg(); }
66
67 ✗ friend Self operator ++ (Self& a, int) { Self r = a; ++a; return r; }
68 ✗ friend Self operator -- (Self& a, int) { Self r = a; --a; return r; }
69 23745420542 friend Self operator + (const Self& a, const Self& b) { return Self(a) += b; }
70 23427826403 friend Self operator - (const Self& a, const Self& b) { return Self(a) -= b; }
71 13721935768 friend Self operator * (const Self& a, const Self& b) { return Self(a) *= b; }
72 13626853 friend Self operator / (const Self& a, const Self& b) { return Self(a) /= b; }
73
74 8201 friend bool operator != (const Self& a, const Self& b) { return !(a == b); }
75
76 ✗ friend Self neg(const Self& a) { return a.neg(); }
77 3375072 friend Self inv(const Self& a) { return a.inv(); }
78 };
79
80 // Storage and arithmetic for numbers mod Self::MOD, as a reduced
81 // representative v in [0, MOD) of unsigned type V.
82 // The type provides static MOD (of type V), reduce (value -> representative),
83 // and *=;
84 // everything else is derived here, valid for any MOD up to V's full range
85 // (sums and differences are tracked mod 2^bits, so no headroom is needed).
86 // Hooks may be overridden in the type's own body (e.g. a faster += / -=).
87 template <typename Self, typename V>
88 struct mod_ops : num_ops<Self> {
89 static_assert(std::unsigned_integral<V>);
90 V v;
91
92 struct is_reduced_tag {};
93
94 4840271106 mod_ops() : v(0) {}
95 42580949 mod_ops(V v_, is_reduced_tag) : v(v_) { assert(v < Self::MOD); }
96 1256965914 template <std::integral I> mod_ops(I x) : v(Self::reduce(x)) {}
97
98 85136675 static Self from_reduced(V v) { return Self(v, is_reduced_tag{}); }
99
100 // A negative value reduces via its nonnegative complement: x = -1 - ~x.
101 830025633 static V reduce(std::signed_integral auto x) {
102 using U = std::make_unsigned_t<decltype(x)>;
103
3/6
None:
✓ Branch 2 → 3 taken 63329 times.
✓ Branch 2 → 4 taken 413967 times.
✗ Branch 5 → 6 not taken.
✓ Branch 5 → 7 taken 6169160 times.
unsigned long wala::mod_ops<wala::mod_goldilocks, unsigned long>::reduce<int>(int):
✗ Branch 2 → 3 not taken.
✗ Branch 2 → 5 not taken.
 1228841791 return x < 0 ? V(Self::MOD - 1 - Self::reduce(U(~x))) : Self::reduce(U(x));
104 }
105
106 204326 explicit operator V() const { return v; }
107 44201176 std::make_signed_t<V> balanced() const {
108
6/6
None:
✓ Branch 2 → 3 taken 11556440 times.
✓ Branch 2 → 4 taken 10521129 times.
✓ Branch 6 → 7 taken 11556440 times.
✓ Branch 6 → 8 taken 10521129 times.
wala::mod_ops<wala::modnum<1000000007>, unsigned int>::balanced() const:
✓ Branch 6 → 7 taken 24971 times.
✓ Branch 6 → 9 taken 21067 times.
 44201176 return std::make_signed_t<V>(Self::MOD-v > v ? v : v - Self::MOD);
109 }
110
111
10/16
✗ Branch 2 → 3 not taken.
✓ Branch 2 → 4 taken 42 times.
✗ Branch 3 → 4 not taken.
✓ Branch 3 → 5 taken 20 times.
✗ Branch 4 → 5 not taken.
✓ Branch 4 → 6 taken 42 times.
✗ Branch 5 → 6 not taken.
✓ Branch 5 → 7 taken 20 times.
✗ Branch 6 → 7 not taken.
✓ Branch 6 → 8 taken 60 times.
✓ Branch 9 → 10 taken 376543 times.
✓ Branch 9 → 11 taken 7894 times.
✓ Branch 13 → 14 taken 3971648 times.
✓ Branch 13 → 15 taken 610 times.
✗ Branch 75 → 76 not taken.
✓ Branch 75 → 77 taken 29 times.
 4561372 friend bool operator == (const Self& a, const Self& b) { return a.v == b.v; }
112
17/17
✓ Branch 11 → 12 taken 16696892 times.
✓ Branch 12 → 13 taken 30644307 times.
✓ Branch 13 → 14 taken 41198909 times.
✓ Branch 14 → 15 taken 27251494 times.
✓ Branch 16 → 17 taken 9088297 times.
✓ Branch 17 → 18 taken 29512950 times.
✓ Branch 18 → 19 taken 43193692 times.
✓ Branch 19 → 20 taken 66934098 times.
✓ Branch 20 → 21 taken 44165059 times.
✓ Branch 22 → 23 taken 644060 times.
✓ Branch 23 → 24 taken 644060 times.
✓ Branch 34 → 35 taken 1919469 times.
✓ Branch 35 → 36 taken 1919469 times.
✓ Branch 36 → 37 taken 1359180 times.
✓ Branch 37 → 38 taken 1359180 times.
✓ Branch 40 → 41 taken 1359630 times.
✓ Branch 41 → 42 taken 1359630 times.
 159722383 friend std::ostream& operator << (std::ostream& out, const Self& n) { return out << n.v; }
113 210008724 friend std::istream& operator >> (std::istream& in, Self& n) { int64_t v_; in >> v_; n = Self(v_); return in; }
114
115 ✗ Self& operator ++ () {
116 ✗ ++v;
117 ✗ if (v == Self::MOD) v = 0;
118 ✗ return self();
119 }
120 ✗ Self& operator -- () {
121 ✗ if (v == 0) v = Self::MOD;
122 ✗ --v;
123 ✗ return self();
124 }
125
16/17
wala::mod_ops<wala::modnum<2013265921>, unsigned int>::operator+=(wala::modnum<2013265921> const&):
✓ Branch 8 → 9 taken 202208976 times.
✓ Branch 8 → 10 taken 290893518 times.
✓ Branch 21 → 22 taken 17321299 times.
✓ Branch 21 → 25 taken 4756214 times.
wala::mod_ops<wala::modnum<998244353>, unsigned int>::operator+=(wala::modnum<998244353> const&):
✓ Branch 8 → 9 taken 3067062451 times.
✓ Branch 8 → 10 taken 4443366985 times.
✓ Branch 11 → 12 taken 10199479 times.
✓ Branch 11 → 13 taken 10195153 times.
✓ Branch 17 → 18 taken 15 times.
✓ Branch 17 → 36 taken 15377153 times.
✓ Branch 31 → 32 taken 15 times.
✓ Branch 60 → 61 taken 116650992 times.
✓ Branch 60 → 62 taken 2417711716 times.
✓ Branch 71 → 72 taken 5116137 times.
✓ Branch 71 → 73 taken 6587158 times.
✗ Branch 85 → 86 not taken.
✓ Branch 85 → 87 taken 11703295 times.
 24292903485 Self& operator += (const Self& o) { v = Self::sub_mod_raw(v, Self::MOD - o.v); return self(); }
126
27/27
✓ Branch 5 → 6 taken 27256182 times.
✓ Branch 5 → 7 taken 19702212 times.
✓ Branch 6 → 7 taken 2 times.
✓ Branch 6 → 8 taken 11 times.
✓ Branch 8 → 9 taken 3156854583 times.
✓ Branch 8 → 10 taken 2686247867 times.
✓ Branch 11 → 12 taken 16255453 times.
✓ Branch 11 → 13 taken 10308309 times.
✓ Branch 16 → 17 taken 2361333077 times.
✓ Branch 16 → 18 taken 2730847960 times.
✓ Branch 19 → 20 taken 2358256427 times.
✓ Branch 19 → 21 taken 2733924610 times.
✓ Branch 22 → 23 taken 2355893560 times.
✓ Branch 22 → 24 taken 2736287477 times.
✓ Branch 25 → 26 taken 2736572239 times.
✓ Branch 25 → 27 taken 2360353732 times.
✓ Branch 25 → 29 taken 17332579 times.
✓ Branch 42 → 43 taken 20396714 times.
✓ Branch 42 → 44 taken 26561680 times.
✓ Branch 45 → 46 taken 25868191 times.
✓ Branch 45 → 47 taken 21090203 times.
✓ Branch 54 → 55 taken 1180798987 times.
✓ Branch 54 → 56 taken 1353563721 times.
✓ Branch 57 → 58 taken 1363580696 times.
✓ Branch 57 → 59 taken 1170782012 times.
✓ Branch 70 → 71 taken 417264 times.
✓ Branch 70 → 72 taken 23046762 times.
 16751668963 Self& operator -= (const Self& o) { v = Self::sub_mod_raw(v, o.v); return self(); }
127 6814046 Self& operator /= (const Self& o) { return self() *= o.inv(); }
128
129 // Returns a - b mod MOD, for b in [0, MOD]; wraparound detects the underflow.
130
48/50
✓ Branch 2 → 3 taken 6011068771 times.
✓ Branch 2 → 4 taken 5255072932 times.
✓ Branch 3 → 4 taken 161158 times.
✓ Branch 3 → 5 taken 1486940 times.
✓ Branch 4 → 5 taken 23436729 times.
✓ Branch 4 → 6 taken 23438267 times.
✓ Branch 5 → 6 taken 7406296477 times.
✓ Branch 5 → 7 taken 6393379509 times.
✓ Branch 6 → 7 taken 77802688 times.
✓ Branch 6 → 8 taken 77424247 times.
✓ Branch 7 → 8 taken 12141324 times.
✓ Branch 7 → 9 taken 12141816 times.
✓ Branch 8 → 9 taken 6420600985 times.
✓ Branch 8 → 10 taken 7416283200 times.
✓ Branch 9 → 10 taken 46874996 times.
✗ Branch 9 → 11 not taken.
✓ Branch 11 → 12 taken 10199479 times.
✓ Branch 11 → 13 taken 10195153 times.
✓ Branch 13 → 14 taken 2364098538 times.
✓ Branch 13 → 15 taken 2728082499 times.
✓ Branch 14 → 15 taken 7454549 times.
✓ Branch 14 → 16 taken 7922619 times.
✓ Branch 16 → 17 taken 2366906321 times.
✓ Branch 16 → 18 taken 2738087637 times.
✓ Branch 19 → 20 taken 2358256427 times.
✓ Branch 19 → 21 taken 2733924610 times.
✓ Branch 22 → 23 taken 2357413680 times.
✓ Branch 22 → 24 taken 2752088656 times.
✓ Branch 25 → 26 taken 2731827305 times.
✓ Branch 25 → 27 taken 2360353732 times.
✓ Branch 26 → 27 taken 443749 times.
✓ Branch 26 → 28 taken 4301185 times.
✓ Branch 28 → 29 taken 10 times.
✓ Branch 28 → 30 taken 5 times.
✓ Branch 34 → 35 taken 8671485 times.
✓ Branch 34 → 36 taken 14792541 times.
✓ Branch 39 → 40 taken 20507933 times.
✓ Branch 39 → 41 taken 26450461 times.
✓ Branch 42 → 43 taken 20396714 times.
✓ Branch 42 → 44 taken 26561680 times.
✓ Branch 51 → 52 taken 1181502678 times.
✓ Branch 51 → 53 taken 1352860030 times.
✗ Branch 54 → 55 not taken.
✓ Branch 54 → 56 taken 60 times.
✓ Branch 67 → 68 taken 8388294 times.
✓ Branch 67 → 69 taken 15075732 times.
✓ Branch 68 → 69 taken 5116137 times.
✓ Branch 68 → 70 taken 6587158 times.
✓ Branch 82 → 83 taken 5116137 times.
✓ Branch 82 → 84 taken 6587158 times.
 18435879850 static V sub_mod_raw(V a, V b) { return a < b ? a - b + Self::MOD : a - b; }
131
132
7/10
wala::mod_ops<wala::mod_goldilocks, unsigned long>::neg() const:
✓ Branch 4 → 5 taken 668 times.
✗ Branch 4 → 8 not taken.
wala::mod_ops<wala::modnum<1000000007>, unsigned int>::neg() const:
✓ Branch 4 → 5 taken 2085 times.
✗ Branch 4 → 8 not taken.
wala::mod_ops<wala::modnum<2013265921>, unsigned int>::neg() const:
✓ Branch 4 → 5 taken 164 times.
✗ Branch 4 → 8 not taken.
wala::mod_ops<wala::modnum<998244353>, unsigned int>::neg() const:
✓ Branch 2 → 3 taken 28171870 times.
✓ Branch 2 → 4 taken 4214777 times.
✓ Branch 4 → 5 taken 2553 times.
✓ Branch 4 → 8 taken 79 times.
 32392196 Self neg() const { return from_reduced(v ? Self::MOD - v : 0); }
133 10188753 Self inv() const { return from_reduced(mod_inv_in_range(v, Self::MOD)); }
134
135 private:
136 8727747 Self& self() { return static_cast<Self&>(*this); }
137 };
138
139
13/15
✓ Branch 2 → 3 taken 1326178 times.
✗ Branch 2 → 5 not taken.
✓ Branch 2 → 7 taken 998 times.
✓ Branch 2 → 10 taken 1325939 times.
✓ Branch 6 → 7 taken 1330720 times.
✓ Branch 6 → 8 taken 76 times.
✓ Branch 6 → 18 taken 98 times.
✓ Branch 8 → 9 taken 146 times.
✗ Branch 11 → 12 not taken.
✓ Branch 11 → 13 taken 1648009 times.
✓ Branch 22 → 23 taken 20 times.
✓ Branch 28 → 29 taken 9 times.
✓ Branch 31 → 32 taken 1919469 times.
✓ Branch 51 → 52 taken 1181502678 times.
✓ Branch 51 → 53 taken 1352860030 times.
 4762119637 template <auto MOD_> struct modnum : mod_ops<modnum<MOD_>, std::make_unsigned_t<decltype(MOD_)>> {
140 using Self = modnum;
141 static_assert(MOD_ > 0, "MOD must be positive");
142 using V = std::make_unsigned_t<decltype(MOD_)>;
143 static constexpr V MOD = V(MOD_);
144
145 using base = mod_ops<modnum, V>;
146
45/54
✗ Branch 2 → 3 not taken.
✓ Branch 2 → 4 taken 6169202 times.
✗ Branch 3 → 4 not taken.
✓ Branch 3 → 5 taken 224429150 times.
✓ Branch 4 → 5 taken 2644881 times.
✓ Branch 4 → 6 taken 42 times.
✓ Branch 5 → 6 taken 28 times.
✓ Branch 5 → 7 taken 6169180 times.
✓ Branch 6 → 7 taken 1325960 times.
✓ Branch 6 → 8 taken 52582799 times.
✗ Branch 7 → 8 not taken.
✓ Branch 7 → 9 taken 3304781 times.
✓ Branch 8 → 9 taken 8 times.
✓ Branch 9 → 10 taken 376563 times.
✓ Branch 9 → 11 taken 7894 times.
✓ Branch 10 → 11 taken 43 times.
✓ Branch 11 → 12 taken 21097472 times.
✓ Branch 11 → 13 taken 23567488 times.
✓ Branch 12 → 13 taken 1696 times.
✓ Branch 12 → 14 taken 644060 times.
✓ Branch 13 → 14 taken 3972963 times.
✓ Branch 13 → 15 taken 610 times.
✓ Branch 14 → 15 taken 7454549 times.
✓ Branch 14 → 16 taken 7922619 times.
✗ Branch 15 → 16 not taken.
✓ Branch 15 → 17 taken 22077513 times.
✗ Branch 16 → 17 not taken.
✓ Branch 16 → 18 taken 590 times.
✗ Branch 17 → 18 not taken.
✓ Branch 17 → 19 taken 10644437 times.
✓ Branch 18 → 19 taken 10618909 times.
✓ Branch 18 → 20 taken 11458604 times.
✓ Branch 21 → 22 taken 20 times.
✓ Branch 25 → 26 taken 79 times.
✓ Branch 25 → 27 taken 15 times.
✓ Branch 26 → 27 taken 1359591 times.
✓ Branch 27 → 28 taken 1359136 times.
✓ Branch 28 → 29 taken 10 times.
✓ Branch 28 → 30 taken 5 times.
✓ Branch 30 → 31 taken 59 times.
✓ Branch 31 → 32 taken 17 times.
✓ Branch 32 → 33 taken 59 times.
✗ Branch 40 → 41 not taken.
✓ Branch 40 → 42 taken 6587187 times.
✗ Branch 60 → 61 not taken.
✓ Branch 60 → 62 taken 29 times.
✓ Branch 65 → 66 taken 6587158 times.
✓ Branch 65 → 67 taken 5116137 times.
✓ Branch 68 → 69 taken 5116137 times.
✓ Branch 68 → 70 taken 6587158 times.
✗ Branch 75 → 76 not taken.
✓ Branch 75 → 77 taken 29 times.
✓ Branch 82 → 83 taken 5116137 times.
✓ Branch 82 → 84 taken 6587158 times.
 833823409 using base::base;
147 using base::v;
148 using base::reduce;
149
150
1/1
✓ Branch 24 → 25 taken 17 times.
 433358636 static V reduce(std::unsigned_integral auto x) { return V(x % MOD); }
151
152 22159332 explicit operator std::make_signed_t<V>() const
153 requires (MOD <= V(std::numeric_limits<std::make_signed_t<V>>::max()))
154 {
155
3/4
✓ Branch 3 → 4 taken 220 times.
✗ Branch 3 → 21 not taken.
✓ Branch 106 → 107 taken 17595 times.
✓ Branch 106 → 108 taken 23068 times.
 22159112 return std::make_signed_t<V>(v);
156 }
157
158 14799095392 Self& operator *= (const Self& o) {
159
36/39
✓ Branch 2 → 3 taken 439584292 times.
✓ Branch 2 → 4 taken 192673200 times.
✓ Branch 5 → 6 taken 2661835292 times.
✓ Branch 5 → 7 taken 3207830920 times.
✗ Branch 6 → 7 not taken.
✓ Branch 6 → 8 taken 10644437 times.
✓ Branch 8 → 9 taken 15530640 times.
✓ Branch 8 → 10 taken 25673124 times.
✓ Branch 11 → 12 taken 16255453 times.
✓ Branch 11 → 13 taken 10308309 times.
✓ Branch 12 → 13 taken 321979 times.
✓ Branch 12 → 14 taken 322072 times.
✓ Branch 13 → 14 taken 1487405 times.
✓ Branch 13 → 15 taken 160526 times.
✓ Branch 15 → 16 taken 494484 times.
✓ Branch 15 → 17 taken 22687400 times.
✗ Branch 17 → 18 not taken.
✓ Branch 17 → 19 taken 10644437 times.
✓ Branch 18 → 19 taken 10618909 times.
✓ Branch 18 → 20 taken 11458604 times.
✓ Branch 25 → 26 taken 2731827305 times.
✓ Branch 25 → 27 taken 2360353732 times.
✓ Branch 28 → 29 taken 446296308 times.
✓ Branch 28 → 30 taken 448311263 times.
✓ Branch 45 → 46 taken 25868191 times.
✓ Branch 45 → 47 taken 21090203 times.
✓ Branch 49 → 50 taken 3293587 times.
✓ Branch 49 → 52 taken 3293600 times.
✓ Branch 53 → 54 taken 60 times.
✓ Branch 54 → 55 taken 1180798987 times.
✓ Branch 54 → 56 taken 1353563721 times.
✓ Branch 57 → 58 taken 1363580696 times.
✓ Branch 57 → 59 taken 1170782012 times.
✓ Branch 71 → 72 taken 5116137 times.
✓ Branch 71 → 73 taken 6587158 times.
✓ Branch 79 → 80 taken 6587158 times.
✓ Branch 79 → 81 taken 5116137 times.
✗ Branch 85 → 86 not taken.
✓ Branch 85 → 87 taken 11703295 times.
 14787067489 if constexpr (sizeof(V) <= 4) v = V(uint64_t(v) * o.v % MOD);
160 else v = V(__uint128_t(v) * o.v % MOD);
161 9165676 return *this;
162 }
163 };
164
165 struct mod_goldilocks : mod_ops<mod_goldilocks, uint64_t> {
166 using Self = mod_goldilocks;
167 static constexpr uint64_t MOD = 0xffffffff00000001ull;
168 static constexpr uint64_t EPS = -MOD;
169 // We have 2^32 is a primitive 6th root of unity.
170 // Note that omega_8 + omega_8^7 == 2^24 - 2^72 == sqrt(2)
171 // We'll pick the root so that 2^24 - 2^72 is our primitive 384th root of unity.
172 static constexpr uint64_t PRIMITIVE_ROOT = 2717;
173
174 using base = mod_ops<mod_goldilocks, uint64_t>;
175
6/7
✓ Branch 4 → 5 taken 96 times.
✗ Branch 7 → 8 not taken.
✓ Branch 7 → 9 taken 48 times.
✓ Branch 8 → 9 taken 11556440 times.
✓ Branch 8 → 10 taken 10521129 times.
✓ Branch 11 → 12 taken 11556440 times.
✓ Branch 11 → 13 taken 10521129 times.
 44298111 using base::base;
176 using base::reduce;
177
3/3
✓ Branch 6 → 7 taken 166 times.
✓ Branch 6 → 8 taken 16 times.
✓ Branch 6 → 18 taken 20 times.
 74224918 mod_goldilocks() = default;
178 ✗ mod_goldilocks(__int128_t a) : base(a < 0 ? uint64_t(MOD - 1 - __uint128_t(~a) % MOD) : uint64_t(__uint128_t(a) % MOD), is_reduced_tag{}) {}
179 ✗ mod_goldilocks(__uint128_t a) : base(uint64_t(a % MOD), is_reduced_tag{}) {}
180
181 // Avoids the division: any uint64_t is within MOD of reduced.
182 85577 static uint64_t reduce(std::unsigned_integral auto x) {
183 static_assert(sizeof(x) <= 8);
184 22195150 uint64_t a = x;
185
1/6
None:
✗ Branch 2 → 3 not taken.
✓ Branch 2 → 4 taken 8240 times.
unsigned long wala::mod_goldilocks::reduce<unsigned int>(unsigned int):
✗ Branch 2 → 3 not taken.
✗ Branch 2 → 4 not taken.
unsigned long wala::mod_goldilocks::reduce<unsigned long>(unsigned long):
✗ Branch 2 → 3 not taken.
✗ Branch 2 → 4 not taken.
 8240 return a >= MOD ? a - MOD : a;
186 }
187
188 // returns a-b, assuming -MOD <= a-b, e.g. b <= MOD
189 1574790070 static uint64_t sub_mod_raw(uint64_t a, uint64_t b) {
190 #if defined(__x86_64__)
191 // TODO: We could try to write this using intrinsics, but GCC sometimes produces the wrong code.
192 1574790070 uint64_t res_wrapped = a;
193 1574790070 uint64_t adjustment = b;
194 1574790070 asm (
195 // AT&T syntax: SRC DST
196 "sub %[y], %[x]\n\t"
197 // Trick from plonky2 implementation:
198 // After sub, flag CF is set iff we underflowed. We want to correct by EPS == 2^32 - 1 iff C is set.
199 // sbb (subtract with borrow) computes DST <- DST - SRC - CF
200 // Thus, we can use the 32-bit form of sbb on a dummy register to load CF ? EPS : 0.
201 // Here, we'll just reuse the original register holding b.
202 "sbb %k[y], %k[y]\n\t"
203 : [x] "+r"(res_wrapped),
204 [y] "+r"(adjustment)
205 :
206 : "cc"
207 );
208 #else
209 uint64_t res_wrapped = a - b;
210 uint64_t adjustment = (res_wrapped > a) ? EPS : 0;
211 #endif
212 740705267 return res_wrapped - adjustment;
213 }
214
215 // Reduce lo + 2^64 * mi + 2^96 * hi, where hi <= MOD
216 834084803 static uint64_t reduce_u160_raw(uint64_t lo, uint32_t mi, uint64_t hi) {
217 // result = lo - hi + EPS * mi
218 // 0 <= lo <= 2^64 - 1 = MOD + EPS - 1
219 // 0 <= EPS * mi <= (2^32 - 1) * EPS = MOD - 1 - EPS
220 // 0 <= hi <= MOD
221 // -MOD <= lo - hi + EPS * mi <= 2*MOD-2
222 // so we do have some leeway
223 834084803 return sub_mod_raw(sub_mod_raw(lo, hi), MOD-(uint64_t(mi)<<32)+mi);
224 }
225
226 834084803 static uint64_t reduce_u128_raw(__uint128_t v) {
227 834084803 uint64_t hi = uint64_t(v >> 64);
228 834084803 uint64_t lo = uint64_t(v);
229 834084803 uint32_t hi_hi = uint32_t(hi >> 32);
230 834084803 uint32_t hi_lo = uint32_t(hi);
231 834084803 return reduce_u160_raw(lo, hi_lo, hi_hi);
232 }
233
234 821735278 Self& operator *= (Self o) {
235 821735278 v = reduce_u128_raw(__uint128_t(v) * __uint128_t(o.v));
236
1/2
✗ Branch 15 → 16 not taken.
✓ Branch 15 → 17 taken 22077513 times.
 821689685 return *this;
237 }
238 };
239
240 5218 template <typename T> T power(T a, long long b) {
241 5218 assert(b >= 0);
242
22/26
wala::mod_goldilocks wala::power<wala::mod_goldilocks>(wala::mod_goldilocks, long long):
✓ Branch 4 → 5 taken 17024 times.
✓ Branch 4 → 6 taken 11545 times.
✗ Branch 6 → 7 not taken.
✗ Branch 6 → 8 not taken.
✓ Branch 7 → 4 taken 28569 times.
✓ Branch 7 → 8 taken 532 times.
✓ Branch 8 → 9 taken 6880 times.
✓ Branch 8 → 10 taken 5914 times.
✗ Branch 9 → 6 not taken.
✗ Branch 9 → 10 not taken.
✓ Branch 11 → 8 taken 12794 times.
✓ Branch 11 → 12 taken 215 times.
wala::modnum<2013265921> wala::power<wala::modnum<2013265921> >(wala::modnum<2013265921>, long long):
✓ Branch 5 → 6 taken 2128 times.
✓ Branch 5 → 7 taken 8885 times.
✓ Branch 8 → 5 taken 11013 times.
✓ Branch 8 → 9 taken 972 times.
✓ Branch 8 → 10 taken 2477 times.
✓ Branch 11 → 8 taken 2917 times.
✓ Branch 11 → 12 taken 110 times.
wala::modnum<998244353> wala::power<wala::modnum<998244353> >(wala::modnum<998244353>, long long):
✓ Branch 5 → 6 taken 21879 times.
✓ Branch 5 → 7 taken 51945 times.
✓ Branch 8 → 5 taken 73824 times.
✓ Branch 8 → 9 taken 4894 times.
✓ Branch 8 → 10 taken 4131 times.
✓ Branch 11 → 8 taken 5414 times.
✓ Branch 11 → 12 taken 218 times.
 140292 T r = 1; while (b) { if (b & 1) r *= a; b >>= 1; a *= a; } return r;
243 }
244
245 template <typename tag> struct dynamic_modnum : mod_ops<dynamic_modnum<tag>, uint32_t> {
246 using Self = dynamic_modnum;
247
248 private:
249 inline static uint32_t MOD_ = 0;
250 inline static uint64_t BARRETT_M = 0;
251
252 public:
253 // Make only the const-reference public, to force the use of set_mod
254 static constexpr uint32_t const& MOD = MOD_;
255
256 using base = mod_ops<dynamic_modnum, uint32_t>;
257 using base::base;
258 using base::v;
259 using base::reduce;
260
261 // Barret reduction taken from KACTL:
262 /**
263 * Author: Simon Lindholm
264 * Date: 2020-05-30
265 * License: CC0
266 * Source: https://en.wikipedia.org/wiki/Barrett_reduction
267 * Description: Compute $a \% b$ about 5 times faster than usual, where $b$ is constant but not known at compile time.
268 * Returns a value congruent to $a \pmod b$ in the range $[0, 2b)$.
269 * Status: proven correct, stress-tested
270 * Measured as having 4 times lower latency, and 8 times higher throughput, see stress-test.
271 * Details:
272 * More precisely, it can be proven that the result equals 0 only if $a = 0$,
273 * and otherwise lies in $[1, (1 + a/2^64) * b)$.
274 */
275 ✗ static void set_mod(int mod) {
276 ✗ assert(mod > 0);
277 ✗ MOD_ = uint32_t(mod);
278 ✗ BARRETT_M = (uint64_t(-1) / MOD);
279 }
280 ✗ static uint32_t barrett_reduce_partial(uint64_t a) {
281 ✗ return uint32_t(a - uint64_t((__uint128_t(BARRETT_M) * a) >> 64) * MOD);
282 }
283 ✗ static uint32_t barrett_reduce(uint64_t a) {
284 ✗ int32_t res = int32_t(barrett_reduce_partial(a) - MOD);
285 ✗ return uint32_t((res < 0) ? res + int32_t(MOD) : res);
286 }
287
288 struct mod_reader {
289 ✗ friend std::istream& operator >> (std::istream& i, mod_reader) {
290 ✗ int mod; i >> mod;
291 ✗ Self::set_mod(mod);
292 ✗ return i;
293 }
294 };
295 ✗ static mod_reader MOD_READER() {
296 ✗ return mod_reader();
297 }
298
299 ✗ static uint32_t reduce(std::unsigned_integral auto x) {
300 static_assert(sizeof(x) <= 8);
301 ✗ return barrett_reduce(x);
302 }
303
304 ✗ explicit operator int() const { return int(v); }
305
306 ✗ Self& operator *= (const Self& o) {
307 ✗ v = barrett_reduce(uint64_t(v) * o.v);
308 ✗ return *this;
309 }
310 };
311
312 template <typename T> struct mod_constraint {
313 T v, mod;
314
315 2127 friend mod_constraint operator & (mod_constraint a, mod_constraint b) {
316
2/2
✓ Branch 2 → 3 taken 1036 times.
✓ Branch 2 → 4 taken 1091 times.
 2127 if (a.mod < b.mod) std::swap(a, b);
317
2/2
✓ Branch 4 → 5 taken 109 times.
✓ Branch 4 → 6 taken 2018 times.
 2127 if (b.mod == 1) return a;
318
319 2018 extended_gcd_result<T> egcd = extended_gcd<T>(a.mod, b.mod);
320 2018 assert(a.v % egcd.gcd == b.v % egcd.gcd);
321
322
1/2
✗ Branch 13 → 14 not taken.
✓ Branch 13 → 15 taken 2018 times.
 2018 T extra = b.v - a.v % b.mod;
323
1/2
✗ Branch 15 → 16 not taken.
✓ Branch 15 → 17 taken 2018 times.
 2018 extra /= egcd.gcd;
324
325 2018 extra *= egcd.coeff_a;
326
2/4
✗ Branch 17 → 18 not taken.
✓ Branch 17 → 19 taken 2018 times.
✗ Branch 19 → 20 not taken.
✓ Branch 19 → 21 taken 2018 times.
 2018 extra %= b.mod / egcd.gcd;
327
2/2
✓ Branch 21 → 22 taken 680 times.
✓ Branch 21 → 23 taken 1338 times.
 2018 extra += (extra < 0) ? b.mod / egcd.gcd : 0;
328
329 return mod_constraint{
330 2018 a.v + extra * a.mod,
331 2018 a.mod * (b.mod / egcd.gcd)
332 2018 };
333 }
334 };
335
336 } // namespace wala
337