5054. 统计

题目大意

给定n,k,求满足以下条件的整数数组a的数量：

1. a中共有k个元素
2. a[i]∈[1,n]
3. ∀i∈[1,k),a[i]≤a[i+1]
4. gcd(a1,a2...ak)=1

答案取模1e9+7

Data Constraint
n≤109,k≤1000

题解

显然可以莫比乌斯反演。设g(x)表示从1−n 中取出 k 个，gcd 是 x 的倍数
的⽅案数。那么

Ans=∑d=1nμ(d)×g(d)

g(x)如何求？令bi=ai−ai−1，那么显然bi要满足b1>0,bi≥0(i>1),∑bi≤⌊nx⌋，那么b的方案显然就是a的方案。
现在的问题就类似将⌊nx⌋个物品分成k+1组，第一组不为空，其余组任意。所以可以先将一个物品强制分给第一组。然后就是Ck⌊nx⌋+k−1。
然后这个反演显然可以在O(n√)内完成，问题是如何处理μ的前缀和。

---

杜教筛

杜教筛预处理μ的前缀和是其一个经典应用。
记S(n)=∑ni=1μ(i)
首先对于μ具有这个性质，∑d|nμ(d)=[n==1]
所以

∑i=1n∑d|iμ(d)=∑i=1n∑d|iμ(id)=1

枚举约数d，转化一下

∑d=1n∑i=1⌊nd⌋μ(i)=1

提出一个S(n)

S(n)+∑d=2nS(⌊nd⌋)=1

然后就能递归求S(n)了。
对于较小的n可以预处理S(n)，还可以按照⌊Maxnn⌋记忆化。

---

SRC

#include<cstdio>#include<cstdlib>#include<cstring>#include<iostream>#include<algorithm>#include<cmath>using namespace std ;#define N 100000 + 10#define M 10000000 + 10typedef long long ll ;const int MO = 1e9 + 7 ;const int UP = 1e7 ;bool flag[M] ;int fac[N] , _fac[N] ;int Pri[M] , mu[M] , Sum[M] ;int tag[N] , f[N] ;int T , n , k , Con , ans ;int gcd( int x , int y ) { return y == 0 ? x : gcd( y , x % y ) ; }int Power( int x , int k ) {    int s = 1 ;    while ( k ) {        if ( k & 1 ) s = (ll)s * x % MO ;        x = (ll)x * x % MO ;        k /= 2 ;    }    return s ;}void Pre() {    fac[0] = _fac[0] = 1 ;    for (int i = 1 ; i <= 100000 ; i ++ ) {        fac[i] = (ll)fac[i-1] * i % MO ;        _fac[i] = Power( fac[i] , MO - 2 ) ;    }    mu[1] = 1 ;    for (int i = 2 ; i <= UP ; i ++ ) {        if ( !flag[i] ) {            Pri[++Pri[0]] = i ;            mu[i] = -1 ;        }        for (int j = 1 ; j <= Pri[0] ; j ++ ) {            if ( 1ll * i * Pri[j] > UP ) break ;            flag[i*Pri[j]] = 1 ;            if ( i % Pri[j] == 0 ) { mu[i*Pri[j]] = 0 ; break ; }            mu[i*Pri[j]] = -mu[i] ;        }    }    for (int i = 1 ; i <= UP ; i ++ ) Sum[i] = (Sum[i-1] + mu[i] + MO) % MO ;}ll C( int n , int m ) {    if ( n <= 100000 && m <= 100000 ) return (ll)fac[n] * _fac[m] % MO * _fac[n-m] % MO ;    ll ret = _fac[k] ;    for (int i = n - k + 1 ; i <= n ; i ++ ) ret = ret * i % MO ;    return ret ;}ll CalcSum( int x ) {    if ( x <= UP ) return Sum[x] ;    int t = x <= Con ? x : n / x + Con ;    if ( tag[t] == T ) return f[t] ;    ll ret = 1 ;    for (int d = 2 ; d <= x ; d ++ ) {        int r = x / (x / d) ;        ret = ((ret - CalcSum(x/d) * (r - d + 1) % MO) % MO + MO) % MO ;        d = r ;    }    tag[t] = T ;    return f[t] = ret ;}void Mobius() {    ll ret = 0 ;    for (ll d = 1 ; d <= n ; d ++ ) {        int r = n / (n / d) ;        CalcSum(r) ;        ret = ((ret + 1ll * (CalcSum(r) - CalcSum(d-1) + MO) % MO * C( n / d + k - 1 , k ) % MO) % MO + MO) % MO ;        C( n / d + k - 1 , k ) ;        d = r ;    }    printf( "%lld\n" , ret ) ;}int main() {    freopen( "count.in" , "r" , stdin ) ;    freopen( "count.out" , "w" , stdout ) ;    memset( tag , -1 , sizeof(tag) ) ;    scanf( "%d" , &T ) ;    Pre() ;    while ( T -- ) {        scanf( "%d%d" , &n , &k ) ;        Con = sqrt(n) ;        Mobius() ;    }    return 0 ;}

以上.