題意是模擬一個循環，一開始有一個空序列，之後每次循環：

1.從1到m中隨機選出一個數字添加進去，每個數字被選的概率相同。

2.檢查這個序列的gcd是否為1，如果為1則停止，若否則重復1操作直至gcd為1為止。

求這個序列的長度期望。

也是花了一晚上學習了一下期望dp。

設dp[i]表示當前gcd為i，到gcd為1時添加的元素個數期望。

然後就是傳統的期望dp模型了：

dp[i]=∑p[i→j]dp[j]+w[i→j]

此處w為1，因為每次是添加1個元素

初始化狀態dp[1]=0，因為當gcd為1的時候已經無法再添加元素

狀態轉移就是枚舉i的因數j，然後計算1到m中有多少個數字x使得gcd(x,i)=j，設個數為tp，另一方面，還要計算有多少個數字y使得gcd(y,i)=i，設個數為z，從而有：

z=m/i（此處除法為向下取整）

dp[i]=z/m*dp[i]+Σ(tp/m*dp[j])+1 （此處的除法為取模意義下的除法，即乘以逆元）

也就是

dp[i]=(Σ(tp/m*dp[j])+1)*m/(m-z) （除法意義同上）

最後，由於起點並未明確確定，此處要手動設定起點，對於每個起點，都有1/m的概率選到，所以答案就是

1+Σdp[i]/m （取模下除法）

至於求tp，就是對x/i這個數字質因數分解之後容斥定理求個數，由於本人手殘這部分寫掛了好幾次，終於也是在千辛萬苦之後才寫對

  1 #include<bits/stdc++.h>
  2
 
 using namespace std;
  3 #define ll long long
  4 const ll mod=1e9+7;
  5 ll q_p(ll a,ll b)
  6 {
  7     ll ans=1;
  8     while(b)
  9     {
 10         if(b&1)
 11         {
 12             ans*=a;
 13             ans%=mod;
 14         }
 15         b>>=1;
 16         a*=a;
 
 17         a%=mod;
 18     }
 19     return ans;
 20 }
 21 ll inv(ll x)
 22 {
 23     return q_p(x,mod-2);
 24 }
 25 
 26 ll ret;
 27 vector<ll>vec;
 28 void dfs(ll idx,ll dep,ll lim,ll num,ll tmp)
 29 {
 30     if(num>100000) return;
 31     if(dep==lim)
 32     {
 33         if(lim%2)
 34             ret+=tmp/num;
 35         else
 36             ret-=tmp/num;
 37         return;
 38     }
 39     if(idx>=vec.size()) return;
 40     dfs(idx+1,dep+1,lim,num*vec[idx],tmp);
 41     dfs(idx+1,dep,lim,num,tmp);
 42 }
 43 
 44 bool vis[100005];
 45 ll calc(ll x,ll k,ll n)
 46 {
 47     ll tmp=n/k;
 48     ll tt=x/k;
 49     for(ll i=2;;i++)
 50     {
 51         while(tt%i==0)
 52         {
 53             if(!vis[i]) vec.push_back(i),vis[i]=1;
 54             tt/=i;
 55         }
 56         if(i>sqrt(tt)) i=tt-1;
 57         if(tt==1) break;
 58     }
 59     ret=0;
 60     for(int i=1;i<=vec.size();i++)
 61         dfs(0,0,i,1,tmp);
 62     for(int i=0;i<vec.size();i++) vis[vec[i]]=0;
 63     vec.clear();
 64     return tmp-ret;
 65 }
 66 
 67 
 68 ll dp[100005];
 69 int main()
 70 {
 71     #ifdef amori
 72     clock_t start = clock();
 73     #endif //amori
 74 
 75     ll m;
 76     cin>>m;
 77     dp[1]=0;
 78     ll invm=inv(m);
 79     for(ll i=2;i<=m;i++)
 80     {
 81         dp[i]=1;
 82         for(ll j=1;j<=sqrt(i);j++)
 83         {
 84             if(i%j==0)
 85             {
 86                 //cout<<i<<" "<<j<<" "<<calc(i,j,m)<<" "<<calc(i,i/j,m)<<endl;
 87                 dp[i]+=dp[j]*invm%mod*calc(i,j,m);
 88                 dp[i]%=mod;
 89                 if(j!=1 && i!=j*j)
 90                 {
 91                     dp[i]+=dp[i/j]*invm%mod*calc(i,i/j,m);
 92                     dp[i]%=mod;
 93                 }
 94             }
 95         }
 96         ll tp=m/i;
 97         dp[i]=dp[i]*m%mod*inv(m-tp);
 98         dp[i]%=mod;
 99     }
100     ll sum=0;
101     for(int i=1;i<=m;i++)
102     {
103         sum+=dp[i];
104         sum%=mod;
105     }
106     cout<<sum*invm%mod+1<<endl;
107 
108     #ifdef amori
109     clock_t end = clock();
110     cout<<"Done in "<<end-start<<"ms"<<endl;
111     #endif // amori
112 }

是不是寫的很爛，寫的很爛就對了

別人構造級數求和一下就過了，本蒟蒻還在搞期望dp，頂不住鴨。

Codeforces 1139D（期望dp）