您的位置:首页 > 其它

ZOJ 3868(容斥原理+快速幂)

2015-04-27 01:55 260 查看
GCD Expectation
Time Limit: 4 Seconds      Memory Limit: 262144 KB

Edward has a set of n integers {a1, a2,...,an}. He randomly picks a nonempty subset {x1, x2,…,xm}
(each nonempty subset has equal probability to be picked), and would like to know the expectation of [gcd(x1, x2,…,xm)]k.
Note that gcd(x1, x2,…,xm) is the greatest common divisor of {x1, x2,…,xm}.

Input

There are multiple test cases. The first line of input contains an integer T indicating the number of test cases. For each test case:
The first line contains two integers n, k (1 ≤ n, k ≤ 106). The second line contains n integers a1, a2,…,an (1
≤ ai ≤ 106).
The sum of values max{ai} for all the test cases does not exceed 2000000.

Output

For each case, if the expectation is E, output a single integer denotes E · (2n - 1) modulo 998244353.

Sample Input

1
5 1
1 2 3 4 5

Sample Output

42

Author: LIN, Xi
Source: The 15th Zhejiang University Programming Contest
Submit    Status

题意:给出一个序列,求这个序列子集gcd的k次方和的期望

容斥原理。。。

#include<iostream>
#include<algorithm>
#include<cstdio>
#include<vector>
#include<cstring>
#include<map>
using namespace std;
typedef long long ll;
const int maxn = 1e6 + 10;
const ll mod = 998244353;
#define rep(i,a,b) for(int i=(a);i<(b);i++)
#define pb push_back
int cnt[maxn],a[maxn];
ll d[maxn];
ll quick_exp(ll a,int n)
{
ll res = 1;
while(n) {
if(n&1)res = res*a % mod;
a = a*a % mod;
n >>= 1;
}
return res;
}
int main()
{
int T;
scanf("%d",&T);
while(T--) {
int n,k;
scanf("%d%d",&n,&k);
int Mgcd = 1;
for(int i = 0; i < n; i++) {
int x;
scanf("%d",&x);
Mgcd = max(Mgcd,x);
a[i] = x;
}
memset(cnt,0,sizeof(cnt[0])*Mgcd+20);
for(int i = 0; i < n; i++)cnt[a[i]]++;

ll ans = 0;
/*d[i] 表示子集的gcd为i的方案数*/
for(int i = Mgcd; i > 0; i--) {
int tot = 0;
for(int j = i; j <= Mgcd; j+= i)
{
tot += cnt[j];
d[i] = (d[i]-d[j])%mod;/*容斥*/
}
/*任选一个i的倍数的非空子集有2^tot-1种可能*/
d[i] = (d[i]+quick_exp(2,tot)-1)%mod;
d[i] = (d[i]+mod)%mod;
ans += d[i]*quick_exp(i,k);
ans %= mod;
}
cout<<ans<<endl;
}
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息