当前位置: 代码迷 >> 综合 >> Codeforces 235 E Number Challenge(莫比乌斯反演)
  详细解决方案

Codeforces 235 E Number Challenge(莫比乌斯反演)

热度:22   发布时间:2023-12-08 10:38:27.0

题目链接:
Codeforces 235 E Number Challenge
题意:

d(i)ii=1aj=1bk=1cd(ijk)a,b,c[1,2000]

分析:
Ans = i=1aj=1bk=1cd(ijk) = i=1a?ai?j=1b?bj?k=1c?ck?gcd(i,j)=gcd(i,k)=gcd(j,k)=1

对于后面的部分反演可得:
Ans=i=1i=a?ai?dmin(b,c)μ(d)d|j,(i,j)=1?bj?d|k,(i,k)=1?ck?

j=dj',k=dk',j=dj',k=dk',
Ans=i=1i=a?ai?dmin(b,c)μ(d)gcd(i,dj)=1?bdj?gcd(i,dk)=1?cdk?
,
因为 gcd(i,dj)=1 ,我们先保证 gcd(i,d)=1, 然后枚举 j:1bd ,保证 gcd(i,j)=1 ,这样就可以使得 gcd(i,dj)=1 ,累加即可。对于 gcd(i,dk)=1 同样处理。
枚举 i d 的时间复杂度是是 a min(b,c) ,但是枚举 jk 的时间复杂度是 i=bi=1b/i i=ci=1c/i ,所以总的时间复杂度大概是 O(a?b?log(b))

#include <cstdio>
#include <cstring>
#include <string>
#include <cmath>
#include <algorithm>
#include <bitset>
using namespace std;
typedef long long ll;
const int MAX_N = 2010;
const ll mod = (ll)1 << 30;int prime_cnt, prime[MAX_N], mu[MAX_N], gcd[MAX_N][MAX_N];
bitset<MAX_N> bs;void GetMu()
{prime_cnt = 0;bs.set();mu[1] = 1;for(int i = 2; i < MAX_N; ++i) {if(bs[i]) {prime[prime_cnt++] = i;mu[i] = -1;}for(int j = 0; j < prime_cnt && i * prime[j] < MAX_N; ++j) {bs[i * prime[j]] = 0;if(i % prime[j]) {mu[ i * prime[j]] = - mu[i];}else {mu[i * prime[j]] = 0;break;}}   }
}inline int GCD(int a, int b)
{if(gcd[a][b]) return gcd[a][b];else return b == 0 ? a : gcd[a][b] = GCD(b, a % b);
}inline void GetGcd()
{for(int i = 1; i < MAX_N; i++){for(int j = i; j < MAX_N; j++) {gcd[i][j] = gcd[j][i] = GCD(i, j);}}
}inline ll work(int n, int x)
{ll res = 0;for(int i = 1; i <= n; ++i) {if(gcd[i][x] == 1) {res = (res + (n / i)) % mod;}}return res;
}inline ll solve(int a, int b, int c)
{ll res = 0;int top = min(b, c);for(int i = 1; i <= a; ++i) {for(int d = 1; d <= top; ++d){if(gcd[i][d] == 1) {ll tmp = 0;tmp = (ll) (a / i) * mu[d] * work(b / d, i) % mod * work(c / d, i) % mod;res = ((res + tmp) % mod + mod) % mod;}}}return res;
}int main()
{GetMu();GetGcd();int a, b, c;while(~scanf("%d%d%d", &a, &b, &c)) {printf("%lld\n", solve(a, b, c));}return 0;
}
  相关解决方案