做题计划

唐。

Luogu P3455 [POI2007] ZAP-Queries

莫比乌斯反演。

令:ab

求:

i=1aj=1b[gcd(i,j)=x]

消掉 x

=i=1axj=1bx[gcd(i,j)=1]

根据莫比乌斯的定义:

=i=1axj=1bxd|gcd(i,j)μ(d)

改成和的形式:

=i=1axj=1bxd=1axμ(d)×[dgcd(i,j)]

移项:

=d=1axμ(d)i=1bxj=1ax×[dgcd(i,j)]

显然 d|i,d|j,又因为满足条件的 i1ax 中出现 ax×d 次,j1bx 中出现 bx×d 次,所以:

=d=1axμ(d)ax×dbx×d

然后可以预处理出 μ 的前缀和,另外的整除分块即可。时间复杂度 O(mlogn)

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 5e4; int pri[N + 100], cnt, mu[N + 100], a, b, x, n, sumu[N + 100]; bool vis[N + 100]; void solve() { cin >> a >> b >> x; if (a > b) { swap(a, b); } a /= x; b /= x; int ans = 0; for (int l = 1, r; l <= a; l = r + 1) { r = min(a / (a / l), b / (b / l)); ans += (a / l) * (b / l) * (sumu[r] - sumu[l - 1]); } cout << ans << endl; } void init() { mu[1] = 1; vis[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } else { mu[i * pri[j]] = mu[i] * mu[pri[j]]; } } } for (int i = 1; i <= N; i++) { sumu[i] = sumu[i - 1] + mu[i]; } } signed main() { init(); cin >> n; while (n--) { solve(); } return 0; }

Luogu P2261 [CQOI2007] 余数求和

求:

i=1nkmodi

根据 mod 的性质:

amodb=ab×ab

得:

i=1nki×ki

整除分块即可,注意一下边界。

如果 kl,那么直接跳出循环,否则 r=min(n,kll)

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; int n, k, ans; signed main() { cin >> n >> k; for (int l = 1, r; l <= n; l = r + 1) { if (k <= l) { break; } else { r = min(n, k / (k / l)); } ans += 1ll * (r - l + 1) * (k / l) * (l + r) / 2ll; } cout << n * k - ans; return 0; }

Luogu P2257 YY 的 GCD

求:

i=1nj=1m[gcd(i,j)prime]

枚举 k

=k=1,kprimeni=1nj=1m[gcd(i,j)=k]

除以 k

=k=1,kprimeni=1nkj=1mk[gcd(i,j)=1]

莫比乌斯定义转换:

=k=1,kprimeni=1nkj=1mkp|gcd(i,j)μ(p)

套路:

=k=1,kprimenp=1nkμ(p)×nkp×mkp

T=kp

原式:

=k=1,kprimenp=1nkμ(p)×nT×mT

枚举 T

=T=1nnT×mTkT,kprimeμ(Tk)

线性筛预处理 + 整除分块即可。

点击查看代码
#include <bits/stdc++.h> using namespace std; const int N = 1e7 + 100; int n, m, mu[N], cnt, pri[N], f[N], t; long long sumu[N]; bool vis[N]; void init() { mu[1] = 1; vis[1] = 1; for (int i = 2; i <= N - 100; ++i) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N - 100; ++j) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } else { mu[i * pri[j]] = mu[i] * mu[pri[j]]; } } } for (int i = 1; i <= cnt; ++i) { for (int j = 1; j * pri[i] <= N - 100; ++j) { f[j * pri[i]] += mu[j]; } } for (int i = 1; i <= N - 100; ++i) { sumu[i] = sumu[i - 1] + f[i]; // cout << sumu[i] << ' '; } } void solve() { cin >> n >> m; long long sum = 0; for (int l = 1, r; l <= min(n, m); l = r + 1) { r = min(n / (n / l), m / (m / l)); sum += (sumu[r] - sumu[l - 1]) * (n / l) * (m / l); } cout << sum << '\n'; } int main() { ios::sync_with_stdio(0); cin.tie(0), cout.tie(0); init(); cin >> t; while (t--) { solve(); } return 0; }

SPOJ SP26017

求:

i=1nj=1mgcd(i,j)

令:nm

F1

枚举 d,除以 d

=d=1ndj=1ndj=1md[gcd(i,j)=1]

莫比乌斯反演:

=d=1ndj=1ndj=1mdpgcd(i,j)μ(p)

枚举 p

=d=1ndp=1ndi=1nd[pi]j=1md[pj]μ(p)

套路:

=d=1ndp=1ndnpdmpdμ(p)

T=pd

原式:

=d=1ndp=1ndnTmTμ(p)

枚举 T

=T=1np|TnTmTμ(p)×Tp

调转顺序:

=T=1nnTmTp|Tμ(p)×Tp

F2

有:

d|nφ(d)=n

带入原式:

=i=1nj=1md|gcd(i,j)φ(d)

套路:

=d=1ni=1nj=1mφ(d)×[di]×[d|j]

接着套路:

=d=1nφ(d)×nd×md

给出 F2 代码:

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 5e4; const int mod = 1e9 + 7; int n, m, n1, m1, T, pri[N + 100], phi[N + 100], cnt, sum[N + 100], p1, p2; bool vis[N + 100]; int num(int n, int m) { if (n > m) { swap(n, m); } int ans = 0; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); // cout << "l = " << l << ' ' << "r = " << r << '\n'; // cout << n << ' ' << m << ' ' << min(n / (n / l), m / (m / l)) << endl; ans = (ans + (phi[r] - phi[l - 1]) * (n / l) % mod * (m / l) % mod) % mod; } return ans; } void solve() { cin >> n >> m >> n1 >> m1; // cout << num(n1, m1) << ' ' << num(n - 1, m1) << ' ' << num(n1, m - 1) << ' ' << num(n - 1, m - 1) << endl; cout << ((num(n1, m1) - num(n - 1, m1) - num(n1, m - 1) + num(n - 1, m - 1)) % mod + mod) % mod << endl; } void init() { phi[1] = 1; vis[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; phi[i] = i - 1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { phi[i * pri[j]] = phi[i] * pri[j]; break; } else { phi[i * pri[j]] = phi[i] * (pri[j] - 1); } } } for (int i = 1; i <= N; i++) { phi[i] += phi[i - 1]; phi[i] %= mod; // sum[i] = sum[i - 1] + phi[i]; } // cout << sum[100] << endl; } signed main() { init(); cin >> T; cin >> p1 >> p2; while (T--) { solve(); } return 0; }

Luogu P6055 [RC-02] GCD

求:

i=1nj=1np=1njq=1nj[gcd(i,j)=1][gcd(p,q)=1]

左右同乘 j

i=1nj=1np=1nq=1n[gcd(i,j)=1][gcd(p,q)=j]

所以 j=gcd(p,q),带入 gcd(i,j)=1,得:

=i=1np=1nq=1n[gcd(i,gcd(p,q))=1]

根据 gcd 的性质(gcd(a,gcd(b,c))=gcd(a,b,c)):

=i=1np=1nq=1n[gcd(i,p,q)=1]

莫比乌斯反演:

=i=1np=1nq=1nd|gcd(i,p,q)μ(d)

套路:

=d=1nnd3μ(d)

杜教筛维护 μ 函数,整除分块即可。

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 5e5; const int mod = 998244353; int cnt, n, num, sum[N + 100], pri[N + 100]; bool vis[N + 100]; unordered_map<int, int> mu; void init() { mu[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j]) { mu[i * pri[j]] = -mu[i]; } else { mu[i * pri[j]] = 0; break; } } } for (int i = 1; i <= N; i++) { sum[i] = sum[i - 1] + mu[i]; sum[i] %= mod; } } int suma(int n) { if (n <= N) { return sum[n]; } if (mu[n]) { return mu[n]; } int ans = 1ll; for (int l = 2, r; l <= n; l = r + 1) { r = n / (n / l); ans -= (r - l + 1) * suma(n / l); } return mu[n] = ans; } signed main() { ios::sync_with_stdio(0); cin.tie(0); cout.tie(0); cin >> n; init(); for (int l = 1, r; l <= n; l = r + 1) { r = n / (n / l); int w = (n / l) * (n / l) % mod * (n / l) % mod; // cout << l << ' ' << r << endl; num = (num + (suma(r) - suma(l - 1) + mod) % mod * w % mod) % mod; } cout << num; return 0; }

Luogu P4449 于神之怒加强版

求:

i=1nj=1mgcd(i,j)k

枚举 d

=d=1ndki=1nj=1m[gcd(i,j)=k]

除以 d

=d=1ndki=1ndj=1md[gcd(i,j)=1]

莫比乌斯反演:

=d=1ndki=1ndj=1mdp=1ndμ(p)[d|i][d|j]

套路:

=d=1ndkp=1ndμ(p)npdmpd

T=pd

=d=1ndkp=1ndμ(p)nTmT

枚举 T

=T=1nnTmTd|Tμ(Td)dk

线性筛即可。

在线性筛中:

g(x)=d|xμ(xd)dk

显然 g(x)=idk  μ 为范德蒙德卷积),为积性函数(若 f(x),g(x) 均为积性函数,那么 h(x)=d|xf(d)g(xd) 也为积性函数)。

有:

g(x)={(Tx1)kd=Tx1(Tx)kd=Tx0other

那么 g(pk+1)=p(k+1)xpkx=px(pkp(k1)x)=pxg(px)

有两种情况:

fi×prij={fi×fprijgcd(i,prij)=1线fi×prijkprij|iprij1g(prijk+1)

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 5e6; const int mod = 1e9 + 7; int n, m, k, pri[N + 100], phi[N + 100], cnt, T; bool vis[N + 100]; int qpow(int a, int b) { if (!b) { return 1; } int tmp = qpow(a, b / 2); if (b & 1) { return tmp % mod * tmp % mod * a % mod; } return tmp % mod * tmp % mod; } void solve() { cin >> n >> m; if (n > m) { swap(n, m); } int ans = 0; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans = (ans + (phi[r] - phi[l - 1] + mod) % mod * (n / l) % mod * (m / l) % mod) % mod; } cout << ans << endl; return ; } void init() { phi[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; // cout << i << ' ' << k << ' ' << qpow(i, k) << endl; phi[i] = (qpow(i, k) - 1 + mod) % mod; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { phi[i * pri[j]] = phi[i] * qpow(pri[j], k) % mod; break; } else { phi[i * pri[j]] = phi[i] * phi[pri[j]] % mod; } } } for (int i = 1; i <= N; i++) { phi[i] += phi[i - 1]; phi[i] %= mod; // sum[i] = sum[i - 1] + phi[i]; } // cout << phi[10] << endl; } signed main() { ios::sync_with_stdio(0); cin.tie(0), cout.tie(0); cin >> T >> k; init(); while (T--) { solve(); } return 0; }

Luogu P3172 [CQOI2015] 选数

求:

a1=LRa2=LRan=LR[gcd(a1,a2an)=k]

除以 k,令 l=Rk,r=Lk

a1=lra2=lran=lr[gcd(a1,a2an)=1]

枚举 p

a1=lra2=lran=lrp|gcd(a1,a2an)μ(p)

套路:

p=1rμ(p)(rpl1p)n

杜教筛预处理,整除分块即可。

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 1e6; const int mod = 1e9 + 7; int cnt, n, num, sum[N + 100], pri[N + 100], n1, k1, l1, h1; bool vis[N + 100]; unordered_map<int, int> mu; void init() { mu[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j]) { mu[i * pri[j]] = -mu[i]; } else { mu[i * pri[j]] = 0; break; } } } for (int i = 1; i <= N; i++) { // cout << i << ' '; sum[i] = sum[i - 1] + mu[i]; sum[i] %= mod; } } int qpow(int a, int b) { if (!b) { return 1; } int tmp = qpow(a, b / 2); if (b & 1) { return tmp % mod * tmp % mod * a % mod; } return tmp % mod * tmp % mod; } int suma(int n) { if (n <= N) { return sum[n]; } if (mu[n]) { return mu[n]; } int ans = 1ll; for (int l = 2, r; l <= n; l = r + 1) { r = n / (n / l); ans -= (r - l + 1) * suma(n / l); } return mu[n] = ans; } signed main() { ios::sync_with_stdio(0); cin.tie(0); cout.tie(0); cin >> n1 >> k1 >> l1 >> h1; l1--; l1 /= k1; h1 /= k1; // cout << l1 << ' ' << h1 << endl; // 2 2 2 4 init(); for (int l = 1, r; l <= h1; l = r + 1) { if (l1 / l) { r = min(h1 / (h1 / l), l1 / (l1 / l)); } else { r = h1 / (h1 / l); } // cout << l << ' ' << r << endl; num = (num + (suma(r) - suma(l - 1) + mod) % mod * qpow((h1 / l - l1 / l + mod) % mod, n1) % mod) % mod; } cout << num; return 0; }

Luogu P6825 「EZEC-4」求和

毒瘤卡常题。

求:

i=1nj=1ngcd(i,j)i+j

枚举 d

d=1ni=1nj=1ndi+j[gcd(i,j)=d]

除以 d

d=1ni=1ndj=1nddd(i+j)[gcd(i,j)=1]

反演:

d=1ni=1ndj=1nddd(i+j)p|gcd(i,j)μ(p)

套路:

d=1np=1ndμ(p)i=1ndj=1nddd(i+j)[p|i][p|j]

不考虑 [p|i][p|j] 了,直接就是 nd 里面有 npdp 的倍数,但是 i,j 同除多除了一个 p,在 dd(i+j) 中补上一个 pdpd(i+j)

=d=1np=1ndμ(p)i=1npdj=1npddpd(i+j)

T=pd

=d=1np=1ndμ(p)i=1nTj=1nTdT(i+j)

拆开:

=d=1np=1ndμ(p)i=1nT(dT)ij=1nT(dT)j

发现 i,j 两个循环一样:

=d=1np=1ndμ(p)(i=1nT(dT)i)2

发现括号里面为等比数列,有:

s(x)={s(n2)(1+qn2)2|xs(n12)(1+qn+12)2x

O(logn) 求和即可。

点击查看代码
#include <bits/stdc++.h> using namespace std; const int N = 1.5e6; int T, n, p, pri[N + 100], mu[N + 100], cnt; bool vis[N + 100]; namespace FastIO { class FastIOBase { protected: #ifdef OPENIOBUF static const int BUFSIZE=1<<16; char buf[BUFSIZE+1]; int buf_p=0; #endif FILE *target; public: #ifdef OPENIOBUF virtual void flush()=0; #endif FastIOBase(FILE *f): target(f){} ~FastIOBase()=default; }; class FastOutput final: public FastIOBase { #ifdef OPENIOBUF public: inline void flush() { fwrite(buf,1,buf_p,target),buf_p=0; } #endif private: inline void __putc(char x) { #ifdef OPENIOBUF if(buf[buf_p++]=x,buf_p==BUFSIZE)flush(); #else putc(x,target); #endif } template<typename T> inline void __write(T x) { static char stk[64],*top;top=stk; if(x<0) return __putc('-'),__write(-x); do *(top++)=x%10,x/=10; while(x); for(;top!=stk;__putc(*(--top)+'0')); } public: FastOutput(FILE *f=stdout): FastIOBase(f){} #ifdef OPENIOBUF ~FastOutput(){ flush(); } #endif inline FastOutput &operator <<(char x) { return __putc(x),*this; } inline FastOutput &operator <<(const char *s) { for(;*s;__putc(*(s++)));return *this; } inline FastOutput &operator <<(const string &s) { return (*this)<<s.c_str(); } template<typename T> inline enable_if_t<is_integral<T>::value,FastOutput&> operator <<(const T &x) { return __write(x),*this; } template<typename ...T> inline void writesp(const T &...x) { initializer_list<int>{(this->operator<<(x),__putc(' '),0)...}; } template<typename ...T> inline void writeln(const T &...x) { initializer_list<int>{(this->operator<<(x),__putc('\n'),0)...}; } template<typename Iter> inline void writesp(Iter begin,Iter end) { while(begin!=end) (*this)<<*(begin++)<<' '; } template<typename Iter> inline void writeln(Iter begin,Iter end) { while(begin!=end) (*this)<<*(begin++)<<'\n'; } }qout; class FastInput final: public FastIOBase { #ifdef OPENIOBUF public: inline void flush() { buf[fread(buf,1,BUFSIZE,target)]=EOF,buf_p=0; } #endif private: bool __eof,__ok; inline char __getc() { if(__eof) return EOF; #ifdef OPENIOBUF if(buf_p==BUFSIZE) flush(); return ~buf[buf_p]?buf[buf_p++]:(__eof=true,EOF); #else char ch=getc(target); return ~ch?ch:(__eof=true,EOF); #endif } public: FastInput(FILE *f=stdin): FastIOBase(f),__eof(false),__ok(true) #ifdef OPENIOBUF { buf_p=BUFSIZE; } #else {} #endif inline bool eof()const { return __eof; } inline char peek() { return __getc(); } explicit inline operator bool()const { return __ok; } inline FastInput &operator >>(char &x) { while(isspace(x=__getc()));return *this; } template<typename T> inline enable_if_t<is_integral<T>::value,FastInput&> operator >>(T &x) { char ch,sym=0; while(isspace(ch=__getc())); if(ch=='-') sym=1,ch=__getc(); if(__eof) return __ok=false,*this; for(x=0;isdigit(ch);x=(x<<1)+(x<<3)+(ch^48),ch=__getc()); return sym?x=-x:x,*this; } inline FastInput &operator >>(char *s) { while(isspace(*s=__getc())); if(__eof) return __ok=false,*this; for(;!isspace(*s) && !__eof;*(++s)=__getc()); return *s='\0',*this; } inline FastInput &operator >>(string &s) { char str_buf[(1<<8)+1],*p=str_buf; char *const buf_end=str_buf+(1<<8); while(isspace(*p=__getc())); if(__eof) return __ok=false,*this; for(s.clear(),p++;;p=str_buf) { for(;p!=buf_end && !isspace(*p=__getc()) && !__eof;p++); *p='\0',s.append(str_buf); if(p!=buf_end) break; } return *this; } template<typename ...T> inline void read(T &...x) { initializer_list<int>{(this->operator>>(x),0)...}; } template<typename Iter> inline void read(Iter begin,Iter end) { while(begin!=end) (*this)>>*(begin++); } }qin; } // namespace FastIO using FastIO::qin,FastIO::qout; inline long long qpow(int a, int b, int p) { if (!b) { return 1; } long long tmp = qpow(a, b / 2, p); if (b & 1) { return tmp % p * tmp % p * a % p; } return tmp % p * tmp % p; } inline void init() { mu[1] = 1; for (int i = 2; i <= N; ++i) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; ++j) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } mu[i * pri[j]] = -mu[i]; } } } inline long sum(int a, int n, int mod) { if (n == 1) { return a; } if (!(n & 1)) { return 1ll * (1ll + qpow(a, n / 2, mod)) * sum(a, n / 2, mod) % mod; } return (1ll * (1ll + qpow(a, n / 2, mod)) * sum(a, n / 2, mod) % mod + qpow(a, n, mod)) % mod; } int main() { // freopen("1.txt", "w", stdout); qin >> T; init(); while (T--) { qin >> n >> p; long long suma = 0; for (int i = 1; i <= n; ++i) { long long ans = qpow(i, i, p); long long t = ans; for (int j = 1; j <= n / i; ++j) { // suma = (suma + mu[j] * sum(n / (i * j), t, p) % p * sum(n / (i * j), t, p) % p) % p; // t = (t * ans) % p; // qout << sum(t, n / (i * j), p) << endl; suma = 1ll * (suma + (1ll * sum(t, n / (i * j), p) % p * sum(t, n / (i * j), p) % p) * (mu[j] + p) % p) % p; t = t * ans % p; } } qout << suma % p << '\n'; } return 0; }

P1829 [国家集训队] Crash的数字表格 / JZPTAB

求:

i=1nj=1mlcm(i,j)

令:nm

有:

i=1mj=1mijgcd(i,j)

枚举 d

d=1ni=1mj=1m[gcd(i,j)=d]ijd

除以 d

d=1ndi=1ndj=1nd[gcd(i,j)=1]ij

莫比乌斯反演:

d=1ndi=1ndj=1ndp|gcd(i,j)μ(p)ij

移相:

d=1ndp=1ndμ(p)i=1ndi[pi]j=1mdj[pj]

枚举 p 的整数数倍:

d=1ndp=1ndμ(p)up=1ndvp=1mduvp2

提出 p2

d=1ndp=1ndp2μ(p)(u=1npdu)(v=1mpdv)

化一下:

d=1ndp=1ndp2μ(p)S(npd)S(mpd)

其中:

S(n)=n(n+1)2

线性筛预处理,两次整除分块即可。

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 1e7; const int mod = 20101009; int pri[N + 100], mu[N + 100], cnt, n, m, sum[N + 100], ans; bool vis[N + 100]; void init() { mu[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } mu[i * pri[j]] = -mu[i]; } } for (int i = 1; i <= N; i++) { sum[i] = (sum[i - 1] + i * i % mod * mu[i] % mod) % mod; } } int sol(int n) { return (1ll * n * (n + 1) / 2ll) % mod; } int suma(int n, int m) { if (n > m) { swap(n, m); } int ans = 0; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans = (ans + (sum[r] - sum[l - 1] + mod) % mod * sol(n / l) % mod * sol(m / l) % mod) % mod; } return ans % mod; } void solve(int n, int m) { if (n > m) { swap(n, m); } int ans = 0; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans = (ans + 1ll * (r - l + 1) * (l + r) / 2ll % mod * suma(n / l, m / l) % mod) % mod; } cout << ans << endl; } signed main() { init(); cin >> n >> m; solve(n, m); return 0; }

Luogu P3327 [SDOI2015] 约数个数和

求:

i=1nj=1md(ij)

其中:

d(x)=p|x1

有:

d(i×j)=x|iy|j[gcd(x,y)=1]

带入:

i=1nj=1mx|iy|j[gcd(x,y)=1]

套路,消掉 i,j

x=1ny=1mnxmy[gcd(x,y)=1]

反演:

x=1ny=1mnxmyp|gcd(x,y)μ(p)

套路:

x=1nnxy=1mmyp=1n[px][py]μ(p)

xipyjp,满足 xpyp

p=1nμ(p)i=1npj=1mpnipmjp

移项:

p=1nμ(p)i=1npnipj=1mpmjp

设:

f(x)=i=1xxi

有:nip=npi=f(nip)mjp 同理:

p=1nμ(p)f(np)f(mp)

预处理 μ 函数和 f 数组,查询的时候整除分块即可。

时间复杂度 O(NlogN)N=5×104

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 5e4 + 100; int pri[N + 100], mu[N + 100], cnt, t, n, m, sum[N + 100], s[N + 100]; bool vis[N + 100]; void init() { mu[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } mu[i * pri[j]] = -mu[i]; } } for (int i = 1; i <= N; i++) { sum[i] = sum[i - 1] + mu[i]; } for (int i = 1; i <= N; i++) { for (int l = 1, r; l <= i; l = r + 1) { r = i / (i / l); s[i] += (r - l + 1) * (i / l); } } } void solve() { cin >> n >> m; if (n > m) { swap(n, m); } int ans = 0; for (int l = 1, r; l <= n; l = r + 1) { // sqrt(n) r = min(n / (n / l), m / (m / l)); ans += (sum[r] - sum[l - 1]) * s[n / l] * s[m / l]; } cout << ans << endl; return ; } signed main() { init(); cin >> t; while (t--) { solve(); } return 0; }

Luogu P1891 疯狂 LCM & SPOJ LCMSUM - LCM Sum

求:

i=1nlcm(i,n)

根据 lcm 的定义:

i=1ningcd(i,n)

n 放在前面:

ni=1nigcd(i,n)

枚举 d

ni=1ndn[gcd(i,n)=d]id

iid

ndni=1nd[gcd(id,n)=d]id

除以 d

ndni=1nd[gcd(i,nd)=1]i

因为 dn,所以 dnd 对称:

ndni=1d[gcd(i,d)=1]i

有:

ndnφ(d)d2

线性筛预处理即可

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 1e6; int t, n, pri[N + 100], phi[N + 100], cnt, f[N + 100]; bool vis[N + 100]; void init() { phi[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; phi[i] = i - 1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { phi[i * pri[j]] = pri[j] * phi[i]; break; } phi[i * pri[j]] = (pri[j] - 1) * phi[i]; } } for (int i = 1; i <= N; i++) { for (int j = 1; i * j <= N; j++) { if (i == 1) { f[i * j] = 1; } else { f[i * j] += phi[i] * i / 2ll; } } } } signed main() { cin >> t; init(); while (t--) { cin >> n; cout << f[n] * n << '\n'; } return 0; }

Luogu P3911 最小公倍数之和

求:

i=1nj=1nlcm(ai,aj)

定义:

ci=j=1n[aj=i]

转化:

i=1nj=1nci×cj×lcm(i,j)

lcm 的定义:

i=1nj=1nci×cj×ijgcd(i,j)

枚举 d

i=1nj=1nd=1n[gcd(i,j)=d]ci×cj×ijd

除以 d,反演:

d=1ni=1ndj=1ndp|gcd(i,j)μ(p)×cid×cjd×i×j×d

枚举 p

d=1np=1ndμ(p)p2i=1npdj=1npdcipd×cjpd×i×j×d

令:

T=pd

枚举 T

T=1nT(i=1nTi×ciT)2p|Tμ(p)p

前面的暴力算,后面预处理即可。

时间复杂度 O(NlogN)

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 5e4 + 100; int n, a[N], c[N], pri[N], mu[N], cnt, f[N]; bool vis[N]; void init() { mu[1] = 1; for (int i = 2; i <= N - 100; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } mu[i * pri[j]] = -mu[i]; } } for (int i = 1; i <= N - 100; i++) { for (int j = 1; i * j <= N - 100; j++) { f[i * j] += i * mu[i]; } } } int qpow(int a, int b) { if (!b) return 1; int tmp = qpow(a, b / 2); if (b & 1) return tmp * tmp * a; return tmp * tmp; } signed main() { cin >> n; init(); for (int i = 1; i <= n; i++) { cin >> a[i]; c[a[i]]++; } int sum = 0; for (int T = 1; T <= N; T++) { int ans = 0; for (int i = 1; i <= N / T; i++) { ans += i * c[i * T]; } sum += T * qpow(ans, 2) * f[T]; } cout << sum; return 0; }

SPOJ PROD1GCD - Product it again

求:

i=1nj=1mgcd(i,j)

有:

d=1ndi=1nj=1m[gcd(i,j)=d]

看指数部分:

i=1nj=1m[gcd(i,j)=d]

反演:

i=1ndj=1mdp|gcd(i,j)μ(p)

枚举 p

p=1ndnpdmpdμ(p)

代回:

d=1ndp=1ndnpdmpdμ(p)

两次整除分块即可。

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 1e7; const int mod = 1e9 + 7; int t, inv[N + 100], pri[N + 100], mu[N + 100], n, m, cnt; bool vis[N + 100]; void init() { mu[1] = 1; inv[0] = inv[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } mu[i * pri[j]] = -mu[i]; } mu[i] = (mu[i] + mu[i - 1]) % mod; inv[i] = (inv[i - 1] * i) % mod; } } int qpow(int a, int b) { if (!b) { return 1; } int tmp = qpow(a, b / 2); if (b & 1) { return tmp % mod * tmp % mod * a % mod; } return tmp % mod * tmp % mod; } int solve(int n, int m) { if (n > m) { swap(n, m); } int ans = 0; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans = (ans + (mu[r] - mu[l - 1] + (mod - 1)) % (mod - 1) * (n / l) % (mod - 1) * (m / l) % (mod - 1)) % (mod - 1); } return ans % (mod - 1); } signed main() { cin >> t; init(); while (t--) { cin >> n >> m; if (n > m) { swap(n, m); } int ans = 1; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans = (ans * qpow(inv[r] * qpow(inv[l - 1], mod - 2) % mod, solve(n / l, m / l) % mod) % mod) % mod; } cout << ans << endl; } return 0; }

Luogu P2260 [清华集训2012] 模积和

求:

i=1nj=1m(nmodi)×(mmodj)(ij)

容斥:

i=1nj=1m(nmodi)×(mmodj)i=1n(nmodi)×(mmodi)

mod 的定义:

i=1nj=1m(nni×i)(mmj×j)i=1n(nni×i)(mmi×i)

移相:

i=1n(nni×i)j=1m(mmj×j)i=1n(nni×i)(mmi×i)

展开:

(n2i=1ni×ni)(m2i=1mi×mi)i=1n(nmmi×mini×ni+i2×mi×ni)

公式:

i=1ni2=n(n+1)(2n+1)6

点击查看代码
#include <bits/stdc++.h> #define int __int128 using namespace std; const int mod = 19940417; int n, m, ans1, ans2, ans3, p, ans4, ans5; int len(int x) { return x * (x + 1) / 2; } int ipow(int x) { return x * (x + 1) * (2 * x + 1) / 6; } inline int read() { int x = 0, f = 1; char ch = getchar(); while (ch < '0' || ch > '9') { if (ch == '-') { f = -1; ch = getchar(); } } while (ch >= '0' && ch <= '9') { x = x * 10 + ch - '0'; ch = getchar(); } return x * f; } void write(int x) { if (x < 0) { putchar('-'); x = -x; } if (x > 9) { write(x / 10); } putchar(x % 10 + '0'); return ; } signed main() { n = read(); m = read(); if (n > m) { swap(n, m); } for (int l = 1, r; l <= n; l = r + 1) { r = n / (n / l); ans1 = (ans1 + (len(r) - len(l - 1) + mod) % mod * (n / l) % mod) % mod; } ans1 = n * n - ans1; for (int l = 1, r; l <= m; l = r + 1) { r = m / (m / l); ans2 = (ans2 + (len(r) - len(l - 1) + mod) % mod * (m / l) % mod) % mod; } ans2 = m * m - ans2; p = n * n * m; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans3 = (ans3 + (len(r) - len(l - 1) + mod) % mod * (n / l) % mod) % mod; } p -= ans3 * m; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans4 = (ans4 + (len(r) - len(l - 1) + mod) % mod * (m / l) % mod) % mod; } p -= ans4 * n; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans5 = (ans5 + (ipow(r) - ipow(l - 1) + mod) % mod * (n / l) % mod * (m / l)) % mod; } p += ans5; write((ans1 * ans2 % mod - p % mod + mod) % mod); return 0; }

UVA Count LCM

求:

i=1nj=1m[lcm(i,j)=ij]

有:

i=1nj=1m[ijgcd(i,j)=ij]

即求:

i=1nj=1m[gcd(i,j)=1]

反演:

i=1nj=1mp|gcd(i,j)μ(p)

套路:

p=1ni=1n[pi]j=1m[pj]μ(p)

套路:

p=1nnpmpμ(p)

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 1e6; int pri[N + 100], mu[N + 100], cnt, t, n, m; bool vis[N + 100]; void init() { mu[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } mu[i * pri[j]] = -mu[i]; } mu[i] += mu[i - 1]; } } void solve() { cin >> n >> m; if (n > m) { swap(n, m); } int ans = 0; for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans += (mu[r] - mu[l - 1]) * (n / l) * (m / l); } cout << ans << '\n'; return ; } signed main() { cin >> t; init(); while (t--) { solve(); } return 0; }

SPOJ VLATTICE - Visible Lattice Points

转化题意:

i=1nj=1nk=1n[gcd(i,j,k)=1]+3i=1nj=1n[gcd(i,j)=1]+3

反演:

i=1nj=1nk=1nd|gcd(i,j,k)μ(d)+3i=1nj=1nd|gcd(i,j)μ(d)+3

套路:

d=1ni=1n[di]j=1n[dj]k=1n[dk]μ(d)+3d=1ni=1n[di]j=1n[dj]μ(d)+3

套路:

d=1n(nd)3μ(d)+3d=1n(nd)2μ(d)+3

合并同类项:

d=1nμ(d)(nd)2(nd+3)+3

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 1e6; int t, n, mu[N + 100], pri[N + 100], cnt; bool vis[N + 100]; void init() { mu[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } mu[i * pri[j]] = -mu[i]; } mu[i] += mu[i - 1]; } } int qpow(int a, int b) { if (!b) { return 1; } int tmp = qpow(a, b / 2); if (b & 1) { return tmp * tmp * a; } return tmp * tmp; } signed main() { // cout << qpow(1, 1) cin >> t; init(); while (t--) { cin >> n; int ans = 0; for (int l = 1, r; l <= n; l = r + 1) { r = n / (n / l); ans += (mu[r] - mu[l - 1]) * qpow(n / l, 2) * (n / l + 3); } cout << (ans + 3) << endl; } return 0; }

SPOJ GCDMAT2 - GCD OF MATRIX (hard)

此题 的升级版。

发现时间复杂度正确,但是常数过大,考虑优化:

  • 根号分治

  • 把除法变成乘这个数的倒数

稍微卡一卡就可以了。

代码太丑,不放了

Luogu P5221 Product

求:

i=1nj=1nlcm(i,j)gcd(i,j)

有:

i=1nj=1nijgcd(i,j)2

拆开:

i=1nj=1nij×i=1nj=1ngcd(i,j)2

反演:

(n!)2n×(d=1ndp=1ndnpdnpdμ(p))2

点击查看代码
#include <bits/stdc++.h> using namespace std; const int N = 1e6; const long long mod = 104857601; int t, pri[N + 100], mu[N + 100], n, cnt; long long fac; bool vis[N + 100]; void init() { mu[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } mu[i * pri[j]] = -mu[i]; } mu[i] = (mu[i] + mu[i - 1]) % mod; } fac = 1ll; for (int i = 1; i <= n; i++) { fac = (long long)fac * i % mod; } } long long qpow(int a, int b) { if (!b) { return 1; } int tmp = qpow(a, b / 2); if (b & 1) { return tmp % mod * tmp % mod * a % mod; } return tmp % mod * tmp % mod; } int solve(int n) { long long ans = 0; for (int l = 1, r; l <= n; l = r + 1) { r = n / (n / l); ans = (ans + (mu[r] - mu[l - 1] + (mod - 1)) % (mod - 1) * (n / l) % (mod - 1) * (n / l) % (mod - 1)) % (mod - 1); } return ans % (mod - 1); } int main() { cin >> n; init(); long long sum = qpow((long long)fac * fac % mod, n); long long ans = 1ll; for (int l = 1, r; l <= n; l = r + 1) { r = n / (n / l); fac = 1ll; for (int i = l; i <= r; i++) { fac = (long long)fac * i % mod; } long long t = qpow((long long)fac * fac % mod, solve(n / l)); ans = (long long)ans * t % mod; } cout << (long long)sum * qpow(ans, mod - 2) % mod << endl; return 0; }

Atcoder [ABC162E] Sum of gcd of Tuples (Hard)

求:

i=1kia1=1ka2=1kan=1k[gcd({an})=i]

令:

p=ki

有:

i=1kia1=1pa2=1pan=1p[gcd({an})=1]

变成 Luogu P3172 [CQOI2015] 选数(见上文):

i=1kit=1pμ(t)(pt)n

代回:

i=1kit=1kiμ(t)(kit)n

时间复杂度 O(klogk)

点击查看代码
#include <bits/stdc++.h> #define int long long using namespace std; const int N = 1e5; const int mod = 1e9 + 7; int cnt, n, k, pri[N + 100], mu[N + 100], sum; bool vis[N + 100]; void init() { mu[1] = 1; for (int i = 2; i <= N; i++) { if (!vis[i]) { pri[++cnt] = i; mu[i] = -1; } for (int j = 1; j <= cnt && i * pri[j] <= N; j++) { vis[i * pri[j]] = 1; if (i % pri[j]) { mu[i * pri[j]] = -mu[i]; } else { mu[i * pri[j]] = 0; break; } } } } int qpow(int a, int b) { if (!b) { return 1; } int tmp = qpow(a, b / 2); if (b & 1) { return tmp % mod * tmp % mod * a % mod; } return tmp % mod * tmp % mod; } signed main() { cin >> n >> k; init(); for (int i = 1; i <= k; i++) { int ans = 0; for (int j = 1; j <= k / i; j++) { ans = (ans + mu[j] * qpow(k / i / j, n) % mod) % mod; } sum = (sum + i * ans % mod) % mod; } cout << sum << endl; return 0; }

__EOF__

本文作者ようこそ!
本文链接https://www.cnblogs.com/ydq1101/p/17789523.html
关于博主:评论和私信会在第一时间回复。或者直接私信我。
版权声明:本博客所有文章除特别声明外,均采用 BY-NC-SA 许可协议。转载请注明出处!
声援博主:如果您觉得文章对您有帮助,可以点击文章右下角推荐一下。您的鼓励是博主的最大动力!
posted @   ydq1101  阅读(95)  评论(0编辑  收藏  举报
相关博文:
阅读排行:
· 微软正式发布.NET 10 Preview 1:开启下一代开发框架新篇章
· 没有源码,如何修改代码逻辑?
· NetPad:一个.NET开源、跨平台的C#编辑器
· PowerShell开发游戏 · 打蜜蜂
· 凌晨三点救火实录:Java内存泄漏的七个神坑,你至少踩过三个!
点击右上角即可分享
微信分享提示