hdu 4746 Mophues 莫比乌斯反演+前缀和优化

时间:2022-01-19 16:23:44

Mophues

题意:给出n, m, p,求有多少对a, b满足gcd(a, b)的素因子个数<=p,(其中1<=a<=n, 1<=b<=m)

有Q组数据;(n, m, P <= 5×105. Q <=5000).

参考:ACdreamers

思路:对于hdu1695 GCD来说,由于只需要求gcd = k的个数,所以我们可以按照不优化莫比乌斯公式直接求,这样求解的时间复杂度为O(n);

若这题我们也不优化,答案累加需要双重循环:

rep1(i,,n){
if(num[i] <= p){
for(int j = i;j <= n;j += i) //优化对象
ans += mu(c[j/i])*(n/j)*(m/j);
}
}

时间复杂度直接为O(n^2),其中n = min(n,m);num[i]表示i的素因子个数;

优化:对于hdu 1695我们是按照公式2即枚举d,然后得到d的倍数j,乘以mu[]的sigma和;这时发现后面的F[j]一直在变,这是因为我们求解的是特定的f[d].但是在这一题中求的不是对于特定的gcd == d的个数,而是一个区间的个数(并且还是gcd的素因子的个数的区间)。我们就需要看对于每一个F[n]的总系数Σmu[]为多少?(这样就可以线性处理了)如果没有素因子的限制,我们就可以直接在预处理出mu[i]之后,枚举i的倍数j,累加即可;现在有系数了,其实也就是只多了一个维度,我们将j的因子的素因子个数加到对应的位置([j][num[i]]),这样再前缀求和就可以弄出素因子个数在p以内了;

原本以为预处理出了每个j的Σmu后即可不再优化了,即每组数据线性处理,但是数据组数过多,导致还是TLE了。这里还需要用到分块的思想;即再一次计算前缀和mu[1...j][k],因为对于d ε [i,n/(n/i)],有n/d = n/i;这样我们不需要枚举i了,直接整块处理即可;

时间复杂度:预处理是O(nlog(n))(也是空间复杂度);之后每组数据时间复杂度为O(sqrt(min(n,m))),分块加速;

#include<iostream>
#include<cstdio>
#include<cstring>
#include<string.h>
#include<algorithm>
#include<vector>
#include<cmath>
#include<stdlib.h>
#include<time.h>
#include<stack>
#include<set>
#include<map>
#include<queue>
using namespace std;
#define rep0(i,l,r) for(int i = (l);i < (r);i++)
#define rep1(i,l,r) for(int i = (l);i <= (r);i++)
#define rep_0(i,r,l) for(int i = (r);i > (l);i--)
#define rep_1(i,r,l) for(int i = (r);i >= (l);i--)
#define MS0(a) memset(a,0,sizeof(a))
#define MS1(a) memset(a,-1,sizeof(a))
#define MSi(a) memset(a,0x3f,sizeof(a))
#define inf 0x3f3f3f3f
#define lson l, m, rt << 1
#define rson m+1, r, rt << 1|1
typedef pair<int,int> PII;
#define A first
#define B second
#define MK make_pair
typedef __int64 ll;
template<typename T>
void read1(T &m)
{
T x=,f=;char ch=getchar();
while(ch<''||ch>''){if(ch=='-')f=-;ch=getchar();}
while(ch>=''&&ch<=''){x=x*+ch-'';ch=getchar();}
m = x*f;
}
template<typename T>
void read2(T &a,T &b){read1(a);read1(b);}
template<typename T>
void read3(T &a,T &b,T &c){read1(a);read1(b);read1(c);}
template<typename T>
void out(T a)
{
if(a>) out(a/);
putchar(a%+'');
}
const int N = 5e5+;
const int M = ;
int num[N],c[N],mu[N][M];
int cal(int n,int p)
{
int cnt = ;
while(n%p == ){
cnt++;
n /= p;
}
return cnt;
}
void getprime()
{
for(int i = ;i < N;i++)if(!num[i]){//素数;
for(int j = i;j < N;j += i){
int t = cal(j,i);
num[j] += t;
if(t > ) c[j] = -;
else if(c[j] >= ) c[j]++;
}
}
}
int mobius(int n)
{
if(n == -) return ;
if(n & ) return -;
return ;
}
void init()
{
getprime();
rep0(i,,N){
for(int j = i; j < N;j += i){
mu[j][num[i]] += mobius(c[j/i]);// 对于每一个j的因数i,最后累加时都有一个mu[j/i]*F[j];
}
}
rep0(i,,N){
rep0(j,,M)
mu[i][j] += mu[i][j-];// 计算在p范围内有多少有效
}
rep0(j,,M){
rep0(i,,N)
mu[i][j] += mu[i-][j];
}
}
int main()
{
init();
int n,m,p,i,j,T;
read1(T);
while(T--){
read3(n,m,p);
if(p >= M){
printf("%I64d\n",1LL*n*m);
continue;
}
ll ans = ;
if(n > m) swap(n,m);
for(i = ;i <= n;i = j + ){//分块处理
j = min(n/(n/i),m/(m/i));//j代表的是值,并不是公约数;  
ans += 1LL*(mu[j][p] - mu[i-][p])*(n/i)*(m/i);
}
printf("%I64d\n",ans);
}
return ;
}