CF1139D Steps to One 題解【莫比烏斯反演】【枚舉】【DP】
反演套 DP 的好題(不用反演貌似也能做
Description
Vivek initially has an empty array \(a\) and some integer constant \(m\).
He performs the following algorithm:
- Select a random integer \(x\) uniformly in range from \(1\) to \(m\) and append it to the end of \(a\).
- Compute the greatest common divisor of integers in \(a\)
.- In case it equals to \(1\), break
- Otherwise, return to step \(1\).
Find the expected length of \(a\). It can be shown that it can be represented as \(\frac PQ\) where \(P\) and \(Q\) are coprime integers and \(Q\neq 0\pmod{10^9+7}\). Print the value of \(P\cdot Q^{-1}\pmod{10^9+7}\).
Input
The first and only line contains a single integer \(m\)
(\(1\le m\le 100000\)).Output
Print a single integer — the expected length of the array \(a\) written as \(P\cdot Q^{-1}\pmod{10^9+7}\).
Examples
input
1
output
1
input
2
output
2
input
4
output
333333338
Note
In the first example, since Vivek can choose only integers from \(1\) to \(1\), he will have \(a=[1]\)
after the first append operation, and after that quit the algorithm. Hence the length of \(a\) is always \(1\), so its expected value is \(1\) as well.In the second example, Vivek each time will append either \(1\) or \(2\), so after finishing the algorithm he will end up having some number of \(2\)‘s (possibly zero), and a single \(1\) in the end. The expected length of the list is \(1?\frac 12+2?\frac 1{2^2}+3?\frac 1{2^3}+\dots =2\).
題意
每一步在 \(1\sim m\) 中任選一個整數,問期望多少步後選出的數的最大公約數是 \(1\)。答案對 \(1\ 000\ 000\ 007\) 取模。
題解
因為每一步不會讓已經選了的元素的 \(\gcd\) 和變大,因此認為是一個除自環外的有向無環圖。對於自環我們很好處理,所以把它看成是一道期望 DP。
令 \(f[i]\) 表示當前的 \(\gcd\) 和為 \(i\),到 \(\gcd\) 和為 \(1\) 的狀態的期望步數。因此把狀態轉移方程寫出來
\[
f[i]=1+\frac{\sum_{j=1}^m f[\gcd(i,j)]}m
\]
這樣的轉移是 \(O(m^2)\) 的。但是我們發現,對於很多 \(j\),\(\gcd(i,j)\) 都是相等的,因此我們把這樣的數整合到一起。
令 \(F(n)\) 表示 \(1\sim m\) 中有多少個數 \(i\) 滿足 \(\gcd(x,i)=n\),其中視 \(x?\) 為常數。
則計算 \(f[i]\) 就轉化為了
\[
f[i]=1+\frac{\sum_{d|i}f[d]\times F(d)}n,x=i
\]
這樣差不多就把枚舉優化到了 \(\log n\) 的 \(d|i?\)。
考慮怎麽計算 \(F(n)\)
\[
F(n)=\sum_{i=1}^m[\gcd(x,i)=n]
\]
令 \(G(n)=\sum_{n|d}F(d)\),則
\[
\begin{aligned}
G(n)&=\sum_{n|d}F(d)\&=\sum_{n|d}\sum_{i=1}^m[\gcd(x,i)=d]\&=\sum_{i=1}^m[n|\gcd(x,i)]\?&=\sum_{i=1}^{\left\lfloor\frac mn\right\rfloor}\left[1|\gcd\left(\frac xn,i\right)\right]
\end{aligned}
\]
實際上這樣是有問題的,因為(在後面)無法保證 \(n|x\),此時 \(G(n)\) 就一定為 \(0\) 了。
我們再多化一步:
\[
\begin{aligned}
G(n)&=\sum_{i=1}^{\left\lfloor\frac mn\right\rfloor}\left[1|\gcd\left(\frac xn,i\right)\right][n|x]\&=\left\lfloor\frac mn\right\rfloor\cdot[n|x]
\end{aligned}
\]
根據 \(G(n)=\sum_{n|d}F(d)\),我們反演到 \(F\),得
\[
\begin{aligned}
F(n)&=\sum_{n|d}\mu\left(\frac dn\right)G(d)\&=\sum_{i=1}^{\left\lfloor\frac mn\right\rfloor}\mu(i)G(ni)\&=\sum_{i=1}^{\left\lfloor\frac mn\right\rfloor}\mu(i)\left\lfloor\frac{m}{ni}\right\rfloor[(ni)|x]
\end{aligned}
\]
我們發現後面的布爾表達式可以當作條件。原本的條件本來就是 \(\to +\infty\) 的,只不過超過了 \(\left\lfloor\frac mn\right\rfloor\) 沒有意義。因此直接把條件換成 \([(ni)|x]\) 即可。又因為 \(n|x\) 在上面的枚舉過程中是成立的,同時可以轉化為 \(\left[i|\frac xn\right]\)。
\[
F(n)=\sum_{i|\frac xn}\mu(i)\left\lfloor\frac{m}{ni}\right\rfloor
\]
這樣的一次枚舉是 \(O\left(d\left(\frac xn\right)\right)\) 的,由於 \(1\sim m\) 的約數個數和均攤是 \(O(\log m)\) 的,其中最多的有 \(128\) 個約數,但是這樣的數肯定不是很多,並且其中很多被枚舉到的數都是質因數,叠代一下並不會造成很大的復雜度。
然後我們需要再把狀態轉移方程稍微轉化一下,把 \(f[i]\) 移到左邊
\[
\begin{aligned}
f[i]&=1+\frac{\sum_{d|i,d<i}f[d]\times F(d)+f[i]\times F(i)}n,x=i\\frac{n-F(i)}{n}\cdot f[i]&=1+\frac{\sum_{d|i,d<i}f[d]\times F(d)}n,x=i\f[i]&=\frac{n+\sum_{d|i,d<i}f[d]\times F(d)}{n-F(i)},x=i
\end{aligned}
\]
就得到了真正的轉移方程。
時間復雜度 \(O(m\log^2 m)\)。
Code:
#include<cstdio>
#include<cstring>
#include<vector>
#define ll long long
#define p 1000000007
using std::vector;
vector<int> v[100100];//約數用 vector 存一下,每次 √m 枚舉不是很穩
ll qpow(ll x,ll y)
{
ll ans=1;
while(y)
{
if(y&1)
ans=ans*x%p;
x=x*x%p;
y>>=1;
}
return ans;
}
bool is[100100];
int pri[100100],mu[100100],cnt=0;
ll f[100100];
int n;
int calc(int x,int y)//1~n 中 gcd(x,i)=y 的數的個數
{
int g=x/y,ans=0;
for(int i=0;i<v[g].size();++i)
ans+=mu[v[g][i]]*(n/v[g][i]/y);
return ans;
}
int main()
{
scanf("%d",&n);
f[1]=1;
mu[1]=1;
for(int i=2;i<=n;++i)
{
if(!is[i])
{
pri[++cnt]=i;
mu[i]=-1;
}
for(int j=1;j<=cnt&&i*pri[j]<=n;++j)
{
is[i*pri[j]]=1;
if(i%pri[j])
mu[i*pri[j]]=-mu[i];
else
{
mu[i*pri[j]]=0;
break;
}
}
}
for(int i=1;i<=n;++i)
for(int j=i;j<=n;j+=i)
v[j].push_back(i);
ll ans=1,inv=qpow(n,p-2);
for(int i=2;i<=n;++i)
{
for(int j=0;j<v[i].size()-1;++j)
f[i]=(f[i]+calc(i,v[i][j])*f[v[i][j]]%p)%p;
f[i]=(f[i]*inv+1)%p;
ll g=n-calc(i,i);
f[i]=f[i]*n%p*qpow(g,p-2)%p;
ans=(ans+f[i])%p;
}
printf("%lld\n",ans*qpow(n,p-2)%p);
return 0;
}
CF1139D Steps to One 題解【莫比烏斯反演】【枚舉】【DP】