HDU4746 Mophues(莫比烏斯反演)
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×10
5. 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) = (1)*B(a,b,d) + (-1)*B(a,b,2*d)+ (-1) *B(a,b,3*d) +(0)*B(a,b,4*d)+(-1)*B(a,b,5*d)+(1)*B(a,b,6*d)...........
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
#include
#include
#include
#include
#include
#include
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 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<> ncase;
while(ncase--){
scanf("%d%d%d",&n,&m,&p);
if(p>=19){
cout<<(ll)n*m< m) swap(n,m);
solve();
}
return 0;
}
