ZOJ 3868 GCD Expectation（容斥+快速幂）

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

题意：给你一个n(n<=1e6)个数的序列a(a<=1e6),从中任意挑一些数构成子集，求所有子集的gcd的k次方的和mod998244353。

输入：

T

n k

a1~an

思路：

设dp[i]表示gcd为i的子集的个数，那么如果能整除i的数共有sum个，gcd为i的倍数的子集就有2^sum-1个。

所以只需要2^sum-1-dp[i*2]-dp[i*3]-...-dp[i*x]（i*x<=max{a1,a2,...,an}，i*(x+1)>max{a1,a2,...,an}）。

因为a<=1e6,因此这个i*x可以直接枚举。然后再用快速幂计算答案即可。

枚举gcd时从大到小（因为上面那个式子从大转移到小的，也算是一个dp吧。。。）

代码：

#include<bits/stdc++.h>
#define ll long long
#define inf 0x3f3f3f3f
#define rep(i,a,b) for(register int i=(a);i<=(b);i++)
#define dep(i,a,b) for(register int i=(a);i>=(b);i--)
using namespace std;
const int maxn=1e6+5;
//const double pi=acos(-1.0);
//const double eps=1e-9;
const ll mo=998244353;
int n,m,k;
int a[maxn];
int sum[maxn];
int ans,tmp,cnt;
int flag;
ll dp[maxn];
ll power(ll a,ll n)
{
    ll sum=1;
    while(n)
    {
        if(n&1) sum=sum*a%mo;
        n>>=1;
        a=a*a%mo;
    }
    return sum;
}
int main()
{
    int T,cas=1;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&k);
        memset(sum,0,sizeof(sum));
        memset(dp,0,sizeof(dp));
        ans=0; flag=0;
        rep(i,1,n)
        {
            scanf("%d",&a[i]);
            sum[a[i]]++;
            flag=max(flag,a[i]);
        }
        ll ans=0;
        dep(i,flag,1)
        {
            ll cnt=sum[i];
            for(int j=i+i;j<=flag;j+=i)
            {
                cnt+=sum[j];
                dp[i]=(dp[i]-dp[j]+mo)%mo;
            }
            dp[i]=(dp[i]+power(2LL,cnt)-1+mo)%mo;
            ans=(ans+power(i,k)*dp[i]%mo)%mo;
        }
        printf("%lld\n",ans);
    }
    return 0;
}