设(f[i])表示当前集合(a)的(gcd)为(i)的状态下,期望加入多少个数之后(gcd)变为(1)
一点也不显然有
[f[i]=1+frac{1}{m}sumlimits_{j=1}^{m}f[gcd(i,j)]
]
不是很 能过,改为枚举(d=gcd(i,j))
[f[i]=1+frac{1}{m}sumlimits_{d|i}f[d]*cnt(d,i)
]
其中(cnt(d,i)=sumlimits_{j=1}^{m}[gcd(i,j)==d])
设(g(i,j)=sumlimits_{k=1}^{i}[k⊥j])
那么(cnt(d,i)=g(lfloorfrac{m}{d} floor,frac{i}{d}))
考虑求所有(ijle m) 对应的所有(g(lfloorfrac{m}{i} floor,j))
设
[F(a,x)=sumlimits_{x|b}g(a,b)=lfloorfrac{a}{x}
floor
]
根据莫反
[g(a,b)=sumlimits_{x|b}mu(x)lfloorfrac{a}{x}
floor
]
那么求(g(a,b))的时间复杂度就是(O(D(b)))
算一下总时间复杂度
[sumlimits_{i=1}^{m}sumlimits_{j=1}^{lfloorfrac{m}{i}
floor}D(j)=sumlimits_{i=1}^{m}sumlimits_{j=1}^{lfloorfrac{m}{i}
floor}lfloorfrac{lfloorfrac{m}{i}
floor}{j}
floor=sumlimits_{i=1}^{m}lfloorfrac{m}{i}
floor logm=O(mlog^2m)
]
回到原式
[f[i]+1+frac{1}{m}sumlimits_{d|i}f[d]*cnt(d,i)
]
发现(d=i)时无法转移,变化一下
[(1-frac{cnt(i,i)}{m})f[i]=1+frac{1}{m}sumlimits_{d|i,d
eq i}f[d]*cnt(d,i)
]
[f[i]=frac{1}{(1-frac{cnt(i,i)}{m})}(1+frac{1}{m}sumlimits_{d|i,d
eq i}f[d]*cnt(d,i))
]
最后答案为
[frac{1}{m}sumlimits_{i=1}^{m}(1+f[i])
]
#include<bits/stdc++.h>
using namespace std;
namespace red{
#define int long long
#define ls(p) (p<<1)
#define rs(p) (p<<1|1)
#define mid ((l+r)>>1)
#define lowbit(i) ((i)&(-i))
inline int read()
{
int x=0;char ch,f=1;
for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
if(ch=='-') f=0,ch=getchar();
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=3e5+10,p=1e9+7,inf=0x3f3f3f3f;
int n,ans;
int inv[N],ret[N];
int prime[N],u[N],tot;
int num[N],wt[N],f[N];
bool vis[N];
vector<int> div[N],cont[N];
inline int fast(int x,int k)
{
int ret=1;
while(k)
{
if(k&1) ret=ret*x%p;
x=x*x%p;
k>>=1;
}
return ret;
}
inline int solve(int n,int m)
{
int ret=0;
for(auto x:div[m]) ret+=u[x]*(n/x);
return ret;
}
inline void oula()
{
vis[0]=vis[u[1]=1]=1;
for(int i=2;i<=n;++i)
{
if(!vis[i]) u[prime[++tot]=i]=-1;
for(int j=1;j<=tot;++j)
{
if(prime[j]*i>n) break;
vis[prime[j]*i]=1;
if(i%prime[j]==0) break;
u[i*prime[j]]=-u[i];
}
}
for(int i=1;i<=n;++i)
for(int j=i;j<=n;j+=i) div[j].push_back(i);
for(int k,r,l=1;l<=n;)
{
k=n/l,r=n/k;
for(int i=1;i<=k;++i) ret[i]=solve(k,i);
for(int i=l;i<=r;++i)
{
for(int j=i;j<=n;j+=i)
{
cont[j].push_back(ret[j/i]);
}
}
l=r+1;
}
}
inline void main()
{
n=read();inv[1]=1;
for(int i=2;i<=n;++i) inv[i]=(p-p/i)*inv[p%i]%p;
oula();
f[1]=1;//cont[i][j]=cnt(i,j);
for(int i=2;i<=n;++i)
{
int sum=div[i].size();
for(int j=0;j<sum;++j)
{
wt[j+1]=inv[n]*cont[i][j]%p,num[j+1]=div[i][j];
}
int qt=fast((1-wt[sum]+p)%p,p-2);
for(int j=1;j<sum;++j) wt[j]=wt[j]*qt%p;
f[i]=qt;
for(int j=1;j<sum;++j)
f[i]=(wt[j]*f[num[j]]+f[i])%p;
}
for(int i=1;i<=n;++i) ans=(ans+inv[n]*f[i])%p;
printf("%lld
",ans);
}
}
signed main()
{
red::main();
return 0;
}