Skip to content

Instantly share code, notes, and snippets.

@key-moon

key-moon/inv.cpp

Last active Jul 8, 2020
Embed
What would you like to do?
#include <bits/stdc++.h>
using ll = long long;
/*
* task: calculate sum_{i=1}^input inv(i)
*
* input = 20000000
*
* MOD_A := 2^29+11 (smallest prime after 2^29(few bits))
* MOD_B := 2^30-35 (largest prime before 2^30(many bits))
* MOD_C := 2^30+3 (smallest prime after 2^30(few bits))
*
* environment: AtCoder Custom Test(2020/07/08)
* result:
* extgcd
* MOD_A: 3767.0ms (3764, 3767, 3767, 3767, 3767(ms))
* MOD_B: 3774.0ms (3769, 3771, 3775, 3776, 3777(ms))
* MOD_C: 3765.3ms (3763, 3765, 3765, 3766, 3771(ms))
* pow
* MOD_A: 1660.0ms (1657, 1659, 1660, 1661, 1661(ms))
* MOD_B: 2553.3ms (2550, 2552, 2554, 2554, 2554(ms))
* MOD_C: 1932.0ms (1929, 1931, 1932, 1933, 1934(ms))
* MOD_C(specialized)
* : 1730.7ms (1729, 1729, 1730, 1733, 1734(ms))
*/
#define USE_EXTGCD
//smallest prime after 2^29
const int MOD_A = (1 << 29) + 11;
//largest prime before 2^30
const int MOD_B = (1 << 30) - 35;
//smallest prime after 2^30
const int MOD_C = (1 << 30) + 3;
template <const int MOD>
class ModTest{
#ifdef USE_EXTGCD
ll getinv(ll a) {
ll div, p = MOD, x1 = 1, y1 = 0, x2 = 0, y2 = 1;
while (true) {
if (p == 1) return x2 + MOD;
div = a / p; x1 -= x2 * div; y1 -= y2 * div; a %= p;
if (a == 1) return x1 + MOD;
div = p / a; x2 -= x1 * div; y2 -= y1 * div; p %= a;
}
}
#else
ll getinv(ll a) {
if (MOD == MOD_C) {
ll tmp = a;
for (int i = 0; i < 30; i++) tmp = (tmp * tmp) % MOD;
return (a * tmp) % MOD;
}
else {
ll k = MOD - 2, res = 1, tmp = a;
while (k) {
if (k & 1) res = (res * tmp) % MOD;
tmp = (tmp * tmp) % MOD;
k >>= 1;
}
return res;
}
}
#endif
public:
ll Test(ll upper) {
//to avoid optimization
ll res = 0;
for (ll i = 1; i <= upper; i++){
res += getinv(i);
}
return res;
}
};
int main() {
ll upper;
std::cin >> upper;
//ll res = ModTest<MOD_A>().Test(upper);
//ll res = ModTest<MOD_B>().Test(upper);
ll res = ModTest<MOD_C>().Test(upper);
std::cout << res << std::endl;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment