两个结论
gcd(fib(a),fib(b))=fib(gcd(a,b))
GCD FIB
gcd(2^a−1,2^b−1) = 2 ^ gcd(a,b)-1
这个我翻了一本初等数论书然后在第六页就找到了。。。。。。
具体来说通过
gcd(2^a−1,2^b−1)
=gcd(2^a - 1 - (2 ^ b - 1) * 2 ^ (a - b) , 2 ^ b - 1)
=gcd(2 ^ (a - b) - 1 , 2 ^ b - 1)
就可以像辗转相减一样证明了。
然后,枚举d,容斥求gcd(a1,a2,a3…ak)=d && sigma(ai) = n的解
AC Code:
#include<bits/stdc++.h>
#define mod 1000000007
#define maxn 2000005
using namespace std;
int fac[maxn],inv[maxn],invf[maxn],n,k;
inline int C(int n,int m){ return 1ll * fac[n] * invf[m] % mod * invf[n-m] % mod; }
int dp[maxn],F[maxn];
inline int cut(int a){ return a>=mod ? a-mod : a; }
inline int cut2(int a){return a>=(mod-1) ? a - (mod-1) : a; }
inline int Pow(int base,int k){ int ret=1;for(;k;k>>=1,base=1ll*base*base%mod) if(k&1) ret=1ll*ret*base%mod; return ret;}
int main()
{
fac[0]=fac[1]=inv[0]=inv[1]=invf[0]=invf[1]=F[1]=1;
for(int i=2;i<maxn;i++)
fac[i] = 1ll * fac[i-1] * i % mod,
inv[i] = 1ll * (mod - mod / i) * inv[mod % i] % mod,
invf[i] = 1ll * invf[i-1] * inv[i] % mod,
F[i] = cut2(F[i-1] + F[i-2]);
int T;
ios::sync_with_stdio(false);
for(cin>>T;T--;)
{
cin>>n>>k;
int ans=0;
for(int i=n;i>=1;i--)
{
dp[i] = 0;
if(n%i==0)
{
dp[i] = C(n/i+k-1,k-1);
for(int j=i+i;j<=n;j+=i)
dp[i] = cut(dp[i] - dp[j] + mod) % mod;
ans = cut(ans + 1ll * (Pow(2,F[i])-1) * dp[i] % mod);
}
}
printf("%lld\n",1ll*ans*Pow(C(n+k-1,k-1),mod-2)%mod);
}
}