CF1139D Step to one_牛客博客

没想到每日一题也会有黑题

期望DP

设 $dp[i]$ 为当前 $gcd$ 为 $i$ 时需要走 $dp[i]$ 步才能到达1,有

$dp[i] = \sum_{j=1}^{i} dp[gcd(i,j)] / m +1\\ = \sum_{d|i}dp[d] \times f(i , d) / m + 1\\ = \sum_{d | i and d \not= i} dp[d] \times f(i , d) / m + dp[i] \times f(i , i) / m + 1 \\ dp[i] = \frac{n + \sum_{d | i and d \not= i} dp[d] \times f(i , d) / m}{n - \sum_{j=1}^{}f(i , j)} +1$

其中 $f(i,d)$ 表示 $[0,i]$ 中有 $f(i,d)$ 个 $x$ 使得 $gcd(x , i) = d$ ,可以容斥处理

时间复杂度 $O(nlog^2n)$

CODE

#include <algorithm>
#include <cctype>
#include <cmath>
#include <complex>
#include <cstdio>
#include <cstring>
#include <deque>
#include <functional>
#include <list>
#include <map>
#include <iomanip>
#include <iostream>
#include <queue>
#include <set>
#include <stack>
#include <string>
#include <vector>
using namespace std;
#define int long long
#define I inline
#define ri register int
#define For(i , x , y) for(ri i = x ; i <= y ; ++ i)
#define Next(i , x) for(ri i = head[x] ; i ; i = e[i].nxt)
#define pb(x) push_back(x)

I int read() {
    int s = 0 , w = 1; char ch = getchar();
    while(ch < 48 || ch > 57) {if(ch == '-') w = -1; ch = getchar();}
    while(ch >= 48 && ch <= 57) s = (s << 1) + (s << 3) + (ch ^ 48) , ch = getchar();
    return s * w;
}

const int N = 1e5 + 5 , P = 1e9 + 7;

int n , f[N];
vector <int> d[N] , cnt[N];

I int qpow(int x , int y) {
    int res = 1;
    for( ; y ; x = (x * x) % P , y >>= 1)
        if(y & 1) res = (res * x) % P ;
    return res;
}

I int inv(int x) {
    return qpow(x , P - 2);
}


signed main() {
    n = read();
    For(i , 1 , n)
        for(int j = i ; j <= n ; j += i)
            d[j].pb(i) , cnt[j].pb(0);
    For(i , 1 , n)
        for(int j = d[i].size() - 1 ; j >= 0 ; -- j) {
            cnt[i][j] = n / d[i][j];
            for(int k = j + 1 ; k < d[i].size() ; ++ k)
                if(d[i][k] % d[i][j] == 0) cnt[i][j] -= cnt[i][k];
        }
    f[1] = 0;
    For(i , 2 , n) {
        int sum = 0 ,x = 0;
        for(int j = 0 ; j < d[i].size() ; ++ j) {
            if(d[i][j] == i) sum = (sum + cnt[i][j]) % P ;
            else x = (x + f[d[i][j]] * cnt[i][j] % P) % P;
        }
        f[i] = (x + n) % P * inv((n - sum + P) % P) % P;
    }
    int ans = 0;
    For(i , 1 , n) ans = (ans + f[i] + 1) % P;
    cout << ans * inv(n) % P << endl;
    return 0;
}
/*
dp[i] = \sum{ dp[gcd(i , j)] } / m + 1
    = 
*/