hdu 4746 Mophues

Mophues

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 327670/327670 K (Java/Others)
Total Submission(s): 290    Accepted Submission(s): 102

Problem Description

As we know, any positive integer C ( C >= 2 ) can be written as the multiply of some prime numbers:
    C = p1×p2× p3× ... × pk
which p1, p2 ... pk are all prime numbers.For example, if C = 24, then:
    24 = 2 × 2 × 2 × 3
    here, p1 = p2 = p3 = 2, p4 = 3, k = 4

Given two integers P and C. if k<=P( k is the number of C's prime factors), we call C a lucky number of P.

Now, XXX needs to count the number of pairs (a, b), which 1<=a<=n , 1<=b<=m, and gcd(a,b) is a lucky number of a given P ( "gcd" means "greatest common divisor").

Please note that we define 1 as lucky number of any non-negative integers because 1 has no prime factor.

 

Input

The first line of input is an integer Q meaning that there are Q test cases.
Then Q lines follow, each line is a test case and each test case contains three non-negative numbers: n, m and P (n, m, P <= 5×10 ⁵. Q <=5000).

 

Output

For each test case, print the number of pairs (a, b), which 1<=a<=n , 1<=b<=m, and gcd(a,b) is a lucky number of P.

 

Sample Input

  
  

   
   2
10 10 0
10 10 1
  
  

 

Sample Output

  
  

   
   63
93
  
  

 

Source

2013 ACM/ICPC Asia Regional Hangzhou Online

 

Recommend

liuyiding

比赛的时候真是太弱了，这题没做出来，后来看了好多关于gcd(x,y)求(x,y)得对数的统计方法，才知道这个玩意该怎么做。

关于题意：定义F(n)为一个数n的素因子的个数，如12=2*2*3，则F(12)=3，特别的F(1)=0。

然后有10^5个询问，问满足F(gcd(i,j))<=P,(1<=x<=n,1<=y<=m)的(x,y)有多少对？

对于不会数论的人这题真的好难，对于我这种学过数论但是个不会莫比乌斯的弱渣来说更是好难。。。

假设我们枚举gcd(i,j)的值d，即gcd(i,j)=d，如果F(d)<=P，那么我们统计gcd(x,y)=d,(1<=x<=n,1<=y<=m)的(x,y)的对数，再把所有累加起来就是答案，当然了，如果我们可以O(1)得到[gcd(x,y)=d,(1<=x<=n,1<=y<=m)的(x,y)的对数]的话这个算法也是可以接受的，但实际上没有O(1)的方法。。。但我们还是要考虑下这个方法

我们考虑怎样求gcd(x,y)=d，这个问题实际上等价于求gcd(x,y)=1,(1<=x<=[n/d],1<=y<=[m/d])的(x,y)的对数。

然后是怎样求gcd(x,y)=1，

我们设f(d)为gcd(x,y)=d的个数，g(d)为gcd(x,y)是d的倍数的个数，这个值很好计算：g(d)=[n/d]*[m/d]

那么显然有g(d)=f(d)+f(2d)+f(3d)+....

然后由莫比乌斯反演，我们知道f(d)=g(d)*u(1)+g(2*d)*u(2)+g(3*d)*u(3)+....

考虑结果ans，ans为所有f(d)且F(d)<=P的和，

我们枚举1<=i<=n，如果i是某个d的倍数且F(d)<=P，那么ans+=g(i)*u(i/d)=[n/i]*[m/i]*u(i/d)。那么这个怎么计算能更快一点？

我们设G(i)为容斥因子：G(i)=sum{u(i/d) | F(d)<=P} 这个值可以nlogn预处理出来，然后我们只需要ans+=G(i)*[n/i]*[m/i]即可

这样的话总的复杂度为O(n*q)还是会T的样子

然后我们注意到[n/i]*[m/i]在一定的范围内是不变的，这个范围是[i,min(n/(n/i),m/(m/i)]，这样我们可以预处理出G(i)的前缀和，然后加快运算（复杂度网上说是sqrt(n)的。。。）

这样总的复杂度是O(q*sqrt(n)+nlog(n))大概这样，然后就可以过了。。。呵呵呵

#include<cstdio>
#include<cstring>
#include<iostream>
#include<vector>
using namespace std;
const int maxn = 500100;
typedef long long ll;
int mu[maxn],sum[maxn],num[maxn];
ll cnt[maxn][19];
bool flag[maxn];
vector<int>prime;
void init(){
    mu[1]=1;
    for(int i=2;i<maxn;i++){
        if(!flag[i]){
            prime.push_back(i);
            mu[i]=-1;
            num[i]=1;
        }
        for(int j=0;j<prime.size()&&i*prime[j]<maxn;j++){
            flag[i*prime[j]]=true;
            num[i*prime[j]]=num[i]+1;
            if(i%prime[j])mu[i*prime[j]]=-mu[i];
            else {mu[i*prime[j]]=0;break;}
        }
    }
    for(int i=1;i<maxn;i++){
        for(int j=i;j<maxn;j+=i){
            cnt[j][num[i]]+=mu[j/i];
        }
    }
    for(int i=0;i<maxn;i++){
        for(int j=1;j<19;j++){
            cnt[i][j]+=cnt[i][j-1];
        }
    }
    for(int i=1;i<maxn;i++){
        for(int j=0;j<19;j++){
            cnt[i][j]+=cnt[i-1][j];
        }
    }
}
int main(){
    init();
    int q;
    scanf("%d",&q);
    while(q--){
        int n,m,k;
        scanf("%d%d%d",&n,&m,&k);
        k=min(k,18);
        ll ans=0;
        if(n>m)swap(n,m);
        for(int i=1,last=i;i<=n;i=last+1){
            last=min(n/(n/i),m/(m/i));
            ans+=(ll)(cnt[last][k]-cnt[i-1][k])*(n/i)*(m/i);
        }
        //printf("%lld\n",ans);
        printf("%I64d\n",ans);
    }
}