hdu5072 互素数对统计

版权声明:本文为博主原创文章,未经博主允许不得转载。 https://blog.csdn.net/u013514182/article/details/42639665

关键词:统计n个数种与x互素/不互素的个数

解决方案:通过对x素因子分解,得到不同素因子p1,p2...pk,问题转化为n个数中不被p1,p2...pk中任何一个数整除(与x互素)/至少被p1,p2...pk中一个数整除(与x不互素)的个数

心得:通过素因子分解将互素问题转化为整除问题(利用素数非整除即互素的性质)

题意:求n个数中两两不互素或两两互素的三元组数目

称满足上述性质的三元组为P三元组

解决非P三元组的计数问题有两种方法。

法一:先求出数组中与每个a[i]互素的个数f[i],这个最后再求。假设f[i]已知,我们用正难则反的思路求P三元组数目。任意一个非P三元组都既存在互素,又存在不互素的数对,因此我们可以对每个a[i],从数组中找出与a[i]互素和不互素的数各一个,这样就构成了一个非P三元组,而每个非P三元组仅对应两个a[i]选择,因此ans=sum{ f[i]*(n-1-f[i])/2,1<=i<=n }

C(n,3)-ans即为最终结果

法二:既然想到了正难则反的思路,不免要想到用容斥原理。设非P三元组的三个数为a1,a2,a3,Uij表示ai和aj不互素,则非P三元组的个数是:|U12∪U23∪U31|=|U12|+|U23|+|U31|-|U12∩U23|-|U23∩U31|-|U31∩U12|.其中|Uij|=sum{ f[i]*(n-2)*3, 1<=i<=n } |Uij∩Ujk|=sum{ C(f[i],2)*3,1<=i<=n },在上述计算中三元组视为有序,因此最后还要除以3!才是最终结果

最后化简后发现与法一结果相同

最后的最后就是计算f[i]了,这真是个经典的结论,也是用容斥原理来计算的。与欧拉公式的计算过程本质相同,只不过欧拉公式结果可以化简,而本题结果无法化简。具体计算过程见代码  

f[i]求解过程:

已知a[i]的不同质因子分别为p1,p2...pk。Api:a[1,2...n]中是pi的倍数的数集 f[i]=|(!Ap1)∩(!Ap2)∩...∩(!Apk)|

注:容斥原理跑得太慢,本体特别容易超时,这题hdu我都刷了1页半了,对的次数不超过5次 Orz...

#include<stdio.h>
#include<iostream>
#include<string.h>
#include<algorithm>
#include<vector>
#define ll long long
#define sf scanf
#define pf printf
#define INF 1<<29
#define maxn 100010
#define mem(a,b) memset(a,b,sizeof(a))
#define lowbit(x) x&(-x)
const ll mol=1000000007;
using namespace std;

int t,n,a[maxn];
int factor[maxn];
int prime[maxn],cnt;
ll ans;
int is_prime[maxn];
int numprime;
int baseprime[maxn];

void make_prime(){
    for(int i=2;i<=1000;i++) is_prime[i]=1;
    numprime=0;
    for(int i=2;i<1000;i++){
        if(is_prime[i]){
            baseprime[numprime++]=i;
            for(int j=i*i;j<1000;j+=i) is_prime[j]=0;
        }
    }
}

void add(int x){
    for(int i=1;i*i<=x;i++)
        if(x%i==0){
            factor[i]++;
            if(x!=i*i) factor[x/i]++;
        }
}

void divided(int x){
    cnt=0;
    for(int i=0;i<numprime;i++){
        if(baseprime[i]*baseprime[i]>x) break;
        if(x%baseprime[i]==0) prime[cnt++]=baseprime[i];
        while(x%baseprime[i]==0) x/=baseprime[i];
    }
    if(x!=1) prime[cnt++]=x;//加入大于sqrt(x)的素因子!
}

ll solve(){
    ll ans=0;
    for(int i=1;i<(1<<cnt);i++){
        int num=0,tmp=1;
        for(int j=0;j<cnt;j++)
            if((1<<j)&i) { tmp*=prime[j]; num++; }
        if(num&1) ans=ans+factor[tmp];
        else ans=ans-factor[tmp];
    }
    return ans;
}

int main(){
    //freopen("a.txt","r",stdin);
    make_prime();
    scanf("%d",&t);
    while(t--){
        scanf("%d",&n);
        ll N=n;
        ans=0;mem(factor,0);
        for(int i=1;i<=n;i++) { scanf("%d",&a[i]); add(a[i]); }
        for(int i=1;i<=n;i++){
            divided(a[i]);
            ll tmp=solve();
            if(tmp) tmp--;//tmp是与a[i]不互素的数的个数,因为包含a[i]本身,所以要减1,而f[i]=N-1-tmp
            ans=ans+(N-1-tmp)*tmp;
        }
        ans=N*(N-1)*(N-2)/6-ans/2;
        printf("%lld\n",ans);
    }
    return 0;
}


拓展:求与x互素的数的个数


const int maxn = 100000+10;
bool not_prime[maxn];
int n,a[maxn],beinum[maxn],f[maxn];
vector<int> prime;
vector<int> factor[maxn],num[maxn];

void make_prime(){
    for(int i=2;i<=maxn-5;i++){
        if(!not_prime[i]) prime.push_back(i);
        for(int j=0;j<=prime.size()&&i*prime[j]<=maxn-5;j++){
            not_prime[i*prime[j]]=1;
            if(!(i%prime[j])) break;
        }
    }
}

void Factor(int x,int id){
    for(int i=0;i<prime.size()&&prime[i]*prime[i]<=x;i++){
        if(x%prime[i]==0){
            factor[id].push_back(prime[i]);
            int cnt=0;
            while(x%prime[i]==0) { x/=prime[i]; cnt++; }
            num[id].push_back(cnt);
        }
    }
    if(x!=1) { factor[id].push_back(x),num[id].push_back(1); }

}

void jilu(int i){
    int cnt=num[i].size(),tmp;
    for(int j=0;j<(1<<cnt);j++){
        tmp=1;
        for(int k=0;k<cnt;k++) if((1<<k)&j) tmp*=factor[i][k];
        beinum[tmp]++;
    }
}

int main(){
    make_prime();
    scanf("%d",&n);
    for(int i=1;i<=n;i++) { scanf("%d",&a[i]); Factor(a[i],i); jilu(i); }
    for(int i=1;i<=n;i++){
        int cnt=num[i].size();
        for(int j=0;j<(1<<cnt);j++){
            int tmp=1,s=0;
            for(int k=0;k<cnt;k++) if((1<<k)&j) tmp*=factor[i][k],s++;
            if(s&1) f[i]-=beinum[tmp];
            else f[i]+=beinum[tmp];
        }
        if(i==1) f[i]--;
    }
    for(int i=1;i<=n;i++) printf("%d\n",f[i]);
}



展开阅读全文

没有更多推荐了,返回首页