src/nt/mod/montgomery64/lib.hpp
| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | #pragma once | ||
| 2 | |||
| 3 | #include "../../../util/alias/num/lib.hpp" | ||
| 4 | |||
| 5 | namespace tifa_libs { | ||
| 6 | |||
| 7 | template <u64 MOD> | ||
| 8 | struct montgomery64 { | ||
| 9 | static CEXP u64 R = [] { | ||
| 10 | u64 iv = MOD * (2 - MOD * MOD); | ||
| 11 | iv *= 2 - MOD * iv, iv *= 2 - MOD * iv, iv *= 2 - MOD * iv; | ||
| 12 | return iv * (2 - MOD * iv); | ||
| 13 | }(); | ||
| 14 | static CEXP u64 R2 = [] { | ||
| 15 | u64 iv = -MOD % MOD; | ||
| 16 | for (u32 i = 0; i != 64; ++i) | ||
| 17 | if ((iv *= 2) >= MOD) iv -= MOD; | ||
| 18 | return iv; | ||
| 19 | }(); | ||
| 20 | static_assert(MOD & 1); | ||
| 21 | static_assert(R * MOD == 1); | ||
| 22 | static_assert((MOD >> 63) == 0); | ||
| 23 | static_assert(MOD != 1); | ||
| 24 | 1439366060993 | static CEXP u64 mulh(u64 x, u64 y) NE { return u64((u128)x * y >> 64); } | |
| 25 | 709003479263 | static CEXP u64 redc_mul(u64 x, u64 y) NE { | |
| 26 | 709003479263 | u64 res = mulh(x, y) - mulh(x * y * R, MOD); | |
| 27 | 709003479263 | return res + (MOD & -(res >> 63)); | |
| 28 | } | ||
| 29 |
4/14tifa_libs::montgomery64<1000000007ul>::norm(long):
✗ Branch 0 not taken.
✓ Branch 1 taken 20366821 times.
tifa_libs::montgomery64<1012924417ul>::norm(long):
✗ Branch 0 not taken.
✗ Branch 1 not taken.
tifa_libs::montgomery64<167772161ul>::norm(long):
✗ Branch 0 not taken.
✗ Branch 1 not taken.
tifa_libs::montgomery64<2013265921ul>::norm(long):
✓ Branch 0 taken 1199 times.
✗ Branch 1 not taken.
tifa_libs::montgomery64<469762049ul>::norm(long):
✗ Branch 0 not taken.
✗ Branch 1 not taken.
tifa_libs::montgomery64<754974721ul>::norm(long):
✗ Branch 0 not taken.
✗ Branch 1 not taken.
tifa_libs::montgomery64<998244353ul>::norm(long):
✓ Branch 0 taken 6852 times.
✓ Branch 1 taken 563788819 times.
|
584173084 | static CEXP u64 norm(i64 x) NE { return (u64)x + (MOD & u64(-(x < 0))); } |
| 30 | }; | ||
| 31 | template <> // dynamic | ||
| 32 | struct montgomery64<0> { | ||
| 33 | u64 MOD, R, R2; | ||
| 34 | CEXP montgomery64() NE = default; | ||
| 35 | CEXPE montgomery64(u64 m) NE { reset(m); } | ||
| 36 | 4390 | CEXP void reset(u64 m) NE { | |
| 37 |
3/6✓ Branch 0 taken 4390 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4390 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 4390 times.
✗ Branch 5 not taken.
|
4390 | assert(!((m & 1) == 0 || m == 1 || m >> 63)), MOD = m; |
| 38 | 4390 | u64 iv = MOD * (2 - MOD * MOD); | |
| 39 | 4390 | iv *= 2 - MOD * iv, iv *= 2 - MOD * iv, iv *= 2 - MOD * iv, R = iv * (2 - MOD * iv), R2 = -MOD % MOD; | |
| 40 |
2/2✓ Branch 0 taken 280960 times.
✓ Branch 1 taken 4390 times.
|
285350 | flt_ (u32, i, 0, 64) |
| 41 |
2/2✓ Branch 0 taken 152232 times.
✓ Branch 1 taken 128728 times.
|
280960 | if ((R2 *= 2) >= MOD) R2 -= MOD; |
| 42 | 4390 | } | |
| 43 | 193893936747 | ND CEXP u64 mul_h(u64 x, u64 y) CNE { return u64((u128)x * y >> 64); } | |
| 44 | 92288962743 | ND CEXP u64 redc_mul(u64 x, u64 y) CNE { | |
| 45 | 92288962743 | cu64 res = mul_h(x, y) - mul_h(x * y * R, MOD); | |
| 46 | 92288962743 | return res + (MOD & -(res >> 63)); | |
| 47 | } | ||
| 48 | 4568 | ND CEXP u64 norm(i64 x) CNE { return u64(x + i64(MOD & u64(-(x < 0)))); } | |
| 49 | }; | ||
| 50 | |||
| 51 | } // namespace tifa_libs | ||
| 52 |