Codeforces Round #548 (Div. 2) D. Steps to One

https://codeforces.com/contest/1139/problem/D

题意：

给你一个空的序列，现在对其进行操作。

1. 从$[1,m]$中随机取得一个数放到序列的最后。
2. 如果整个序列的$gcd=1$退出。
3. 否则进行第$1$步。

求最后这个序列的期望长度。

首先根据期望的定义：
可以得到 $$len=\sum _{i=1}^{+\infty }\frac{i}{m^i}\cdot R(i)$$
$R(i)$表示当前序列长度为$i$时$gcd=1$并且$i-1$时的$gcd̸=1$的方案数。

假定当前序列长度为$i$：
不考虑$i-1$时的影响时：
令$f(i,n)$表示当前整个序列的$gcd=n$的方案数。
令$g(i,n)$表示当前整个序列的$gcd$是$n$的倍数的方案数。
显然：$$g(i,n)=\sum _{n|d}f(i,d)$$$$g(i,n)=\lfloor \frac{m}{n}{\rfloor }^i$$
莫比乌斯反演可以得到：
$$f(i,n)=\sum _{n|d}g(i,d)\cdot u(\frac{d}{n})$$
只需要知道$f(i,1)$的值
$$f(i,1)=\sum _{d=1}^m\lfloor \frac{m}{d}{\rfloor }^i\cdot u(d)$$
这时需要计算$i-1$的影响：

由上述定义：$$R(i)=f(i,1)-m\cdot f(i-1,1)$$
注意：$R(1)=1$
当$i>1$时满足上述式子。

带入：
$$len=\frac{1}{m}+\sum _{i=2}^{+\infty }\frac{i}{m^i}\sum _{d=1}^m(\lfloor \frac{m}{d}{\rfloor }^i\cdot u(d)-\lfloor \frac{m}{d}{\rfloor }^{i-1}\cdot u(d)\cdot m)$$
化简：
$$len=\frac{1}{m}+\sum _{i=2}^{+\infty }\frac{i}{m^i}\sum _{d=1}^{m}u(d)\cdot \lfloor \frac{m}{d}{\rfloor }^{i-1}\cdot (\lfloor \frac{m}{d}\rfloor -m)$$
交换求和顺序
$$len=\frac{1}{m}+\sum _{d=1}^m(\lfloor \frac{m}{d}\rfloor -m)\cdot u(d)\sum _{i=2}^{+\infty }\frac{i}{m^i}\cdot \lfloor \frac{m}{d}{\rfloor }^{i-1}$$

令 $$q=\frac{1}{m}\cdot \lfloor \frac{m}{d}\rfloor$$
令$$H(n)=\sum _{i=2}{q}^{i-1}\cdot i$$
对$H(n)$ 错位相消：
得到：
$$(1-q)\cdot H(n)=q+\frac{q\cdot (1-q^n)}{1-q}$$
$n=+\infty$ $$(1-q)\cdot H(n)=q+\frac{q}{1-q}$$
$$H(n)=\frac{q\cdot (2-q)}{(1-q{)}^2}$$

$$len=\frac{1}{m}+\sum _{d=1}^m(\lfloor \frac{m}{d}\rfloor -m)\cdot u(d)\cdot \frac{q\cdot (2-q)}{m\cdot (1-q{)}^2}$$

这样就可以在$O(mlogm)$时间内解决。

#include<bits/stdc++.h>
using namespace std;
#define inf 1e18
const int maxn = 4e5+100;
const int mod = 1e9+7;
typedef long long ll;

ll powmod(ll x,ll n) {
    ll ret=1;
    for(;n>0;n>>=1,x=x*x%mod) if(n&1) ret=ret*x%mod;
    return ret;
}
int u[maxn],isp[maxn],tab[maxn];
void init(int n) {
    u[1] = 1;int top =0 ;
    for(int i=2;i<=n;++i) {
        if(!isp[i])tab[top++] = i,u[i]=-1;
        for(int j=0;j<top&&tab[j]*i<=n;++j){
            isp[i*tab[j]] = 1;
            if(i%tab[j]==0) {
                u[i*tab[j]] = 0;
                break;
            }
            u[i*tab[j]] = -u[i];
        }
    }
}
ll cal(ll n,ll d) {
    ll q = n/d*powmod(n,mod-2)%mod;
    ll b = powmod(1-q,mod-2);
    ll ans = q*(2-q)%mod*b%mod*b%mod;
    return (ans + mod)%mod;
}
ll slove(ll n) {
    ll ans = powmod(n,mod-2);
    for(int i=2;i<=n;++i) {
        ans = ans + 1ll*u[i]*((n/i-n)*powmod(n,mod-2)%mod*cal(n,i)%mod);
        ans %= mod;
    }
    return (ans + mod)%mod;
}
int main()
{
    ios::sync_with_stdio(0), cin.tie(0);
    int n,k,cnt=0;
    init(2e5);
    cnt = 0;
    cin>>n;
    cout<<slove(n)<<'\n';
    return 0;
}