hdu多校第六场题解(>=100人)

B - bookshelf
呜呜呜，几个特性不知道，结果式子都没推出来。。
题解链接搓这里
一个放x本书的层美观函数$b(x)=2^{F_x}-1$
整个书架的美观值为各层美观值的最大公约数，考察任意两项:
$gcd(b(x),b(y))$
$=gcd(2^{F_x}-1,2^{F_y}-1)$
$=2^{gcd(F_x,F_y)}-1$
$=2^{F_{gcd(x,y)}}-1$
以上用了gcd的两个特性。
第一个是斐波那契数列的特性，即$gcd(F_x,F_y)=F_{gcd(x,y)}$
第二个是$gcd(c^a+1,c^b+1)$也有特性，即$=c^{gcd(a,b)}+1$
ok,式子推出来之后我们就可以枚举约数，然后对约数容斥即可。
我们可以计算约数至少为g的方案数:
1,每一层都是g的倍数，所以我们只要分配系数即可。
2，系数的和为n/g,共k个系数，每个非负，方案数是个经典组合数$C(n/g+k-1,k-1)$
然后容斥就行啦。

#include<bits/stdc++.h>
using namespace std;
const int mod=1e9+7;
const int phimod=1e9+6;
const int maxn=2e6+5;

long long fmod(long long a)
{
    return a%phimod+((a>=phimod)?phimod:0);
}

//***************************************************
//返回d=gcd(a,b)；和对应于等式ax+by=d中的x，y
long long extgcd(long long a,long long b,long long &x,long long &y)
{
    if(a==0&&b==0)return -1;//无最大公约数
    if(b==0){x=1;y=0;return a;}
    long long d=extgcd(b,a%b,y,x);
    y-=a/b*x;
    return d;//返回gcd(a,b)
}
//****************求逆元******************************
//ax=1(mod n)
long long mod_reverse(long long a,long long n)
{
    long long x,y;
    long long d=extgcd(a,n,x,y);
    if(d==1)return (x%n+n)%n;
    return -1;
}

long long f[maxn];
long long tab1[maxn],tab2[maxn];
void init()
{
    f[1]=1;
    for(int i=2;i<maxn;i++)
        f[i]=fmod(f[i-1]+f[i-2]);
    /*for(int i=1;i<=10;i++)
        cout<<i<<' '<<f[i]<<endl;
    */
    tab1[0]=1;
    for(int i=1;i<maxn;i++)tab1[i]=tab1[i-1]*i%mod;
    tab2[maxn-1]=mod_reverse(tab1[maxn-1],mod);
    for(int i=maxn-2;i>=0;i--)
        tab2[i]=tab2[i+1]*(i+1)%mod;
}

long long quick_mul(long long a,long long b)
{
    long long res=1;
    while(b)
    {
        if(b%2)res=res*a%mod;
        b/=2;
        a=a*a%mod;
    }
    return res;
}
vector<int>V;
int n,k;

long long C(long long a,long long b)
{
    return tab1[a]*tab2[b]%mod*tab2[a-b]%mod;
}
long long cnt[maxn];
void solve()
{
    int e=(int)sqrt(0.0+n);
    for(int i=1;i<=e;i++)if(n%i==0)
    {
        if(i*i==n)V.push_back(i);
        else V.push_back(i),V.push_back(n/i);
    }
    sort(V.begin(),V.end());
    for(int i=0;i<V.size();i++)
    {
        int now=V[i];
        cnt[i]=C(n/now+k-1,k-1);
    }
    long long sum=cnt[0];
    for(int i=V.size()-1;i>=0;i--)
        for(int j=i+1;j<V.size();j++)if(V[j]%V[i]==0)
        cnt[i]=((cnt[i]-cnt[j])%mod+mod)%mod;
    //cout<<quick_mul(8,6)<<endl;
    long long res=0;
    for(int i=0;i<V.size();i++)
    {
        long long tmp=cnt[i];
        tmp=(tmp*(quick_mul(2,f[V[i]])%mod-1)%mod+mod)%mod;
        res=(res+tmp)%mod;
    }
    printf("%lld\n",res*mod_reverse(sum,mod)%mod);
}

int main()
{
    init();
    int T;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d %d",&n,&k);
        V.clear();
        solve();
    }
    return 0;

}