码迷,mamicode.com
首页 > 其他好文 > 详细

HDU4746 Mophues(莫比乌斯反演)

时间:2014-08-27 20:35:28      阅读:374      评论:0      收藏:0      [点我收藏+]

标签:des   style   color   os   java   io   strong   for   ar   

Mophues

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


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×105. 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
 
题意: 5000组样例。 问你[1,n] 和 [1,m]中有多少对数的GCD的素因子个数小于p。

思路:
首先考虑一个相对简单的版本: [1,a] 和 [1,b] 有多少对的数  满足GCD <= d
首先定义两个函数:A(a,b,d) 表示 GCD(a,b) = d的对数,B(a,b,d)表示GCD(a,b) 是d的倍数的对数 易得 B(a,b,d) = (a/d)*(b/d) 根据容斥原理:
A(a,b,d) =   a*b + (-1)*B(a,b,2)+ (-1) *B(a,b,3) +(0)*B(a,b,4)+(-1)*B(a,b,5)+(1)*B(a,b,6)...........
B(a,b,i) 前面的系数正是莫比乌斯函数的值。
那么公式可以写成:
A(a,b,1) =  u(1)*B(a,b,1) + u(2)*B(a,b,2) + u(3) *B(a,b,3) + u(4)*B(a,b,4) + u(5)*B(a,b,5) + u(6)*B(a,b,6)...........
A(a,b,2) =  u(1)*B(a,b,2) + u(2)*B(a,b,4) + u(3) *B(a,b,6) + u(4)*B(a,b,8) + u(5)*B(a,b,10) + u(6)*B(a,b,12)...........
A(a,b,3) =  u(1)*B(a,b,3) + u(2)*B(a,b,6) + u(3) *B(a,b,9) + u(4)*B(a,b,12) + u(5)*B(a,b,15) + u(6)*B(a,b,18)...........
A(a,b,4) =  u(1)*B(a,b,4) + u(2)*B(a,b,8) + u(3) *B(a,b,12) + u(4)*B(a,b,16) + u(5)*B(a,b,20) + u(6)*B(a,b,24)...........
答案就是
A(a,b,1)+A(a,b,2)+A(a,b,3)+......A(a,b,d) =   u(1)*B(a,b,1)+(u(1)+u(2))*B(a,b,2) + ....... (u(1)+u(2)+u(3)+u(6))*B(a,b,6)........
可见A(a,b,d) 前的系数为  sigma(u(i)) (i为d的约数) =  C(a,b,d)

然后,这一题还有一个限制条件,就是要使素因子的个数小于等于p,那么我们定义这个函数D(a,b,d,p) 表示B(a,b,d) 前的系数,那么我们只要从C(a,b,d)中选出一些满足条件的系数即可。 用一个数组F[d][cnt] (cnt为素因子个数)记录,数组表示的是d的因子的素因子个数为cnt的影响因子大小。先计算完单个,再计算前缀和(接下来有用)。
接着,我们发现对于某个d,会满足B(a,b,d) = (B,a,b,d+x),而且  这个 x = min(a/(a/d),b/(b/d)) ,那么整个式子的计算会呈现块状,因此计算这个区间的时候可以用前缀和。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <vector>
#include <string>
#include <algorithm>
#include <queue>
using namespace std;
typedef long long ll;
const int maxn = 5e5+10;
int mobi[maxn];
int preMobi[maxn][20];
int pricnt[maxn];
bool isPrime[maxn];
vector<int> prime;

void getMobi(){
    memset(isPrime,1,sizeof isPrime);
    memset(pricnt,0,sizeof pricnt);
    mobi[1] = 1;
    for(int i = 2; i < maxn; i++){
        if(isPrime[i]){
            mobi[i] = -1;
            pricnt[i] = 1;
            prime.push_back(i);
        }
        for(int j = 0; j < prime.size() && i*prime[j] < maxn; j++){
            pricnt[i*prime[j]] = pricnt[i]+1;
            isPrime[i*prime[j]] = false;
            if(i%prime[j]==0){
                mobi[i*prime[j]] = 0;
                break;
            }else{
                mobi[i*prime[j]] = -mobi[i];
            }
        }
    }
}
void getpreMobi(){
    memset(preMobi,0,sizeof preMobi);
    for(int i = 1; i < maxn; i++){
        for(int j = i; j < maxn; j += i){
            preMobi[j][pricnt[i]] += mobi[j/i];
        }
    }
    for(int i = 1; i < maxn; i++){
        for(int j = 0; j < 20; j++){
            preMobi[i][j] += preMobi[i-1][j];
        }
    }
    for(int i = 0; i < maxn; i++){
        for(int j = 1; j < 20; j++){
            preMobi[i][j] += preMobi[i][j-1];
        }
    }
}
int n,m,p;
ll ans;
void solve(){
    ans = 0;
    for(int i = 1; i <= n; i++){
        int ed = min(n/(n/i),m/(m/i));
        ans += ll(preMobi[ed][p]-preMobi[i-1][p])*(n/i)*(m/i);
        i = ed;
    }
    cout<<ans<<endl;

}
void init(){
    getMobi();
    getpreMobi();

}
int main(){
    init();
    int ncase;
    cin >> ncase;
    while(ncase--){
        scanf("%d%d%d",&n,&m,&p);
        if(p>=19){
            cout<<(ll)n*m<<endl;
            continue;
        }
        if(n > m) swap(n,m);
        solve();
    }
    return 0;
}




 
 

HDU4746 Mophues(莫比乌斯反演)

标签:des   style   color   os   java   io   strong   for   ar   

原文地址:http://blog.csdn.net/mowayao/article/details/38875021

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!