Last active
June 29, 2018 12:26
-
-
Save GoBigorGoHome/c54cf8c6299be3222972fa5475e262a4 to your computer and use it in GitHub Desktop.
hihoCoder #1759 取数字
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#include <bits/stdc++.h> | |
using namespace std; | |
#define pb push_back | |
#define eb emplace_back | |
#define all(x) x.begin(), x.end() | |
#define debug(x) cerr << #x <<": " << (x) << endl | |
//在每个函数的入口处执行一次,出口处执行一次。然后就可以快速得知是哪个地方段错误了 | |
#define DEBUG printf("Passing [%s] in LINE %d\n",__FUNCTION__,__LINE__) | |
#ifdef LOCAL | |
#define see(x) cout << #x << ": " << (x) << endl | |
#endif | |
#ifndef LOCAL | |
#define see(x) | |
#endif | |
#define RED "\033[31m" | |
#define RESET "\033[0m" | |
#define alert(x) cerr << RED << x << RESET << endl | |
#ifndef LOCAL | |
#define see(x) | |
#endif | |
#define Rep(n) for(int _ = 0; _ != (n); ++_) | |
#define rep(i, a, b) for(int i = (a); i <= (b); ++i) | |
#define Rng(i, n) for(int i = 0; i != (n); ++i) | |
#define rng(i, a, b) for(int i = (a); i != (b); ++i) | |
#define RNG(i, a) for(auto &i: (a)) | |
namespace std { | |
template<class T> | |
T begin(std::pair<T, T> p) | |
{ | |
return p.first; | |
} | |
template<class T> | |
T end(std::pair<T, T> p) | |
{ | |
return p.second; | |
} | |
} | |
#if __cplusplus < 201402L | |
template<class Iterator> | |
std::reverse_iterator<Iterator> make_reverse_iterator(Iterator it) | |
{ | |
return std::reverse_iterator<Iterator>(it); | |
} | |
#endif | |
template<class Range> | |
std::pair<std::reverse_iterator<decltype(begin(std::declval<Range>()))>, std::reverse_iterator<decltype(begin(std::declval<Range>()))>> make_reverse_range(Range &&r) | |
{ | |
return std::make_pair(make_reverse_iterator(::begin(r)), make_reverse_iterator(::end(r))); | |
} | |
#define RRNG(x, cont) for (auto &x: make_reverse_range(cont)) | |
template<class T> int sign(const T &a) { return a == 0 ? 0 : a > 0 ? 1 : -1; } | |
template<class T> inline T min(T a, T b, T c){return min(min(a, b), c);} | |
template<class T> inline T max(T a, T b, T c){return max(max(a, b), c);} | |
template<class T> void Min(T &a, const T &b){ a = min(a, b); } | |
template<class T> void Max(T &a, const T &b){ a = max(a, b); } | |
template<typename T> void println(const T &t) { cout << t << '\n'; } | |
template<typename T, typename ...Args> void println(const T &t, const Args &...rest) { cout << t << ' '; println(rest...); } | |
template<typename T> void print(const T &t) { cout << t << ' '; } | |
template<typename T, typename ...Args> void print(const T &t, const Args &...rest) { cout << t; print(rest...); } | |
// this overload is chosen when there's only one argument | |
template<class T> void scan(T &t) { cin >> t; } | |
template<class T, class ...Args> void scan(T &a, Args &...rest) { cin >> a; scan(rest...); } | |
int cas; | |
const double pi = acos(-1); | |
int mod = 1e9 + 7; | |
template<class T> | |
void add_mod(T &a, const T &b) { | |
a += b; | |
if (a >= mod) a -= mod; | |
} | |
template<class T> | |
void sub_mod(T &a, const T &b){ | |
a -= b; | |
if (a < 0) a += mod; | |
} | |
using ll = long long; | |
using ull = unsigned long long; | |
using vec = vector<ll>; | |
using mat = vector<vec>; | |
using pii = pair<int, int>; | |
using pdd = pair<double, double>; | |
using pip = pair<int, pii>; | |
using szt = size_t; | |
using vi = vector<int>; | |
using vl = vector<ll>; | |
using vb = vector<bool>; | |
using vpii = vector<pii>; | |
mat operator*(const mat &a, const mat &b) { | |
mat c(a.size(), vec(b[0].size())); | |
for (int i = 0; i < a.size(); i++) { | |
for (int j = 0; j < a[0].size(); j++) { | |
if (a[i][j]) { // optimization for sparse matrix | |
for (int k = 0; k < b[0].size(); k++) { | |
add_mod(c[i][k], a[i][j] * b[j][k] % mod); | |
} | |
} | |
} | |
} | |
return c; | |
} | |
vec operator*(const mat &a, const vec &b) { | |
vec c(a.size()); | |
for (int i = 0; i < a.size(); i++) { | |
for (int j = 0; j < a[0].size(); j++) { | |
add_mod(c[i], a[i][j] * b[j] % mod); | |
} | |
} | |
return c; | |
} | |
mat pow(mat a, ull n) { | |
mat res(a.size(), vec(a[0].size())); | |
for (int i = 0; i < a.size(); i++) { | |
res[i][i] = 1; | |
} | |
while (n) { | |
if (n & 1) { | |
res = res * a; | |
} | |
a = a * a; | |
n >>= 1; | |
} | |
return res; | |
} | |
ll POW(ll x, ll n) { | |
ll res = 1; | |
for (; n; n /= 2, x *= x, x %= mod) { | |
if (n & 1) { | |
res *= x; | |
res %= mod; | |
} | |
} | |
return res; | |
} | |
ll inv(ll x) { | |
return POW(x, mod - 2); | |
} | |
// 2D rotation | |
void rotate(double &x, double &y, double theta) { | |
double tx = cos(theta) * x - sin(theta) * y; | |
double ty = sin(theta) * x + cos(theta) * y; | |
x = tx, y = ty; | |
} | |
namespace bit { | |
const int BIT_N = 1e5 + 5; | |
int bit[BIT_N]; | |
int sum(int x) { | |
int res = 0; | |
while (x) { | |
res += bit[x]; | |
x -= x & -x; | |
} | |
return res; | |
} | |
int sum(int l, int r) { | |
if (l > r) return 0; | |
return sum(r) - sum(l - 1); | |
} | |
void add(int x, int v, int n) { | |
while (x <= n) { | |
bit[x] += v; | |
x += x & -x; | |
} | |
} | |
} | |
namespace util{ | |
int len(ll x){return snprintf(nullptr, 0, "%lld", x);} | |
vi get_d(ll x){ | |
vi res; | |
while(x) { | |
res.pb(x%10); | |
x /= 10; | |
} | |
reverse(all(res)); | |
return res; | |
} | |
} | |
// #include <ext/pb_ds/priority_queue.hpp> | |
// typedef __gnu_pbds :: priority_queue<pip, less<pip>, __gnu_pbds::thin_heap_tag > Heap; | |
// Heap h; | |
// Heap::point_iterator pos[N][N]; | |
const ll LINF = LLONG_MAX/10; | |
const int INF = INT_MAX/10; | |
const int M = 3000 + 5; | |
const int N = 2e5+5; | |
int main() { | |
// Single Cut of Failure taut me | |
// std::fixed 使得所有实数(默认)输出6位小数,即使实际小数位数多于6位。 | |
cout << setprecision(10) << std::fixed; | |
ios::sync_with_stdio(false); | |
cin.tie(nullptr); | |
#ifdef LOCAL | |
freopen("main.in", "r", stdin); | |
// freopen("main.out", "w", stdout); | |
#endif | |
int n; | |
scan(n, mod); | |
int m = sqrt(n); | |
vi minp(m+1); | |
vi p; | |
p.pb(1); | |
vl f(m+1); | |
f[1] = 1; | |
vi cnt(m+1); // 最小质因子的幂次 | |
vl s1(m+1), s2(m+1), t1(m+1), t2(m+1); | |
vi ub(m+1); // ub[j] 大于j的最小质数在p中的下标 | |
// 1152921504606846976 | |
vl inv(40); | |
inv[1] = 1; // 预处理逆元 | |
for (int i = 2; i < 40; i++) { | |
inv[i]=(mod-mod/i)*inv[mod%i]%mod; | |
} | |
const ll inv2 = (mod+1)/2; | |
rng(i, 2, m+1) { | |
if(!minp[i]) { | |
minp[i] = i; | |
for (int j = i - 1; j >= p.back(); j--) { | |
ub[j] = p.size(); | |
} | |
p.pb(i); | |
cnt[i] = 1; | |
f[i] = inv2; | |
} | |
for(int j=1; j < p.size(); j++) { | |
int x= p[j]; | |
if(x > minp[i]) break; | |
ll t = 1LL * x * i; | |
if(t > m) break; | |
minp[t] = x; | |
if(x < minp[i]) { | |
cnt[t] = 1; | |
f[t] = f[x] * f[i] % mod; | |
} | |
else{ | |
cnt[t] = cnt[i] + 1; | |
ll _t = POW(x, cnt[t]); | |
if(_t == t){ | |
f[t] = inv[cnt[t]+1]; | |
} | |
else{ | |
f[t] = f[t/_t] * f[_t] % mod; | |
} | |
} | |
} | |
} | |
for (int i = m; i >= p.back(); i--) { | |
ub[i] = p.size(); | |
} | |
int t = n / (m+1); | |
for (int i = 1; i <= m; i++) { | |
s1[i] = 1; | |
} | |
// 各种优化 | |
vi v(m+1); | |
for (int i = 1; i <= t; i++) { | |
v[i] = n/i; | |
} | |
for (int i = 1; i <= t; i++) { | |
s2[i] = 1; | |
} | |
auto get_s = [inv2, n, m, &p, &s1, &s2, &ub](int i, ll j)->ll { | |
if (i == p.size() || j < p[i]) return 1; | |
// j >= p[i] | |
if (j < p[i] * p[i]) { | |
int pos = j <= m ? ub[j] : p.size(); | |
return ((pos - i) * inv2 + 1) % mod; | |
} | |
return j <= m ? s1[j] : s2[n / j]; | |
}; | |
for(int i=p.size()-1; i>=1; i--) { | |
int p2 = p[i] * p[i]; | |
for (int k = 1; k <= t; ++k) { | |
int j = v[k]; | |
if (j >= p2) { | |
ll x = p[i]; | |
int c = 1; | |
s2[k] = get_s(i + 1, j); | |
while (j >= x) { //这一部分的复杂度怎么算? | |
s2[k] += get_s(i + 1, j / x) * inv[c + 1] % mod; | |
s2[k] %= mod; | |
x *= p[i]; | |
++c; | |
} | |
} else break; | |
} | |
for (int j = m; j >= p2; j--) { | |
ll k = p[i]; | |
int c = 1; | |
s1[j] = get_s(i + 1, j); | |
while (j >= k) { //这一部分的复杂度怎么算? | |
s1[j] += get_s(i + 1, j / k) * inv[c + 1] % mod; | |
s1[j] %= mod; | |
k *= p[i]; | |
++c; | |
} | |
} | |
} | |
for (int i = 1; i <= t; i++) { | |
t2[i] = v[i]; | |
} | |
for (int i = 1; i <= m; i++) { | |
t1[i] = i; | |
} | |
vi pos1(m+1), pos2(m+1); | |
auto get_t = [inv2, n, m, &pos1, & pos2, &p, &t1, &t2](int i, ll j)->ll{ | |
if (i == 0) return j; | |
if (j <= p[i]) return 1; | |
if (j >= p[i]*p[i]) return j <= m ? t1[j] : t2[n/j]; | |
return j <= m ? t1[j] - (i - pos1[j]) : t2[n/j] - (i - pos2[n/j]); | |
}; | |
for(int i = 1; i < p.size(); i++){ | |
int p2 = p[i] * p[i]; | |
rng(k, 1, t+1) { | |
int j = v[k]; | |
if(j >= p2) { | |
t2[k] = get_t(i-1, j) - get_t(i-1, j/p[i]); | |
pos2[k] = i; | |
} | |
else break; | |
} | |
for (int j = m; j >= p2; j--) { | |
t1[j] = get_t(i-1, j) - get_t(i-1, j/p[i]); | |
pos1[j] = i; | |
} | |
} | |
ll ans = get_s(1, n); | |
for (int i = 1; i <= m; i++) { | |
add_mod(ans, f[i] * (get_t(p.size()-1, n/i) - 1) % mod * inv2 % mod); // 注意:是 p.size()-1 不是 p.size() | |
} | |
println(ans); | |
return 0; | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
实现洲阁筛的过程中,遇到几个问题: