[Codeforces 1139D] Steps to One

[题目链接]

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

[算法]

        考虑dp

        设fi表示现在gcd为i ， 期望多少次gcd变为1

        显然 , fi = (1 / m) * sigma{ fgcd(i , j) } + 1

        直接转移是O(N ^ 2logN)的 ， 显然不能通过

        考虑在转移时枚举gcd ， 显然gcd只可能是i的约数 ， 可以在dp前O(NlogN)预处理每个数的约数

        于是问题转化为求一个形如 : [1 , m]中有多少个数与i的gcd为j的问题

        这等价于求 : [1 , m / j]中有多少个数与(i / j)的gcd为1

        容斥原理计算即可

        时间复杂度 : O(NlogN)( 有较大的常数 )

[代码]

        

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
const int P = 1e9 + 7;
const int MAXP = 1e5 + 10;

#define rint register int

int m , tot;
int f[MAXP] , prime[MAXP] , dp[MAXP];
vector< int > d[MAXP];

template <typename T> inline void chkmax(T &x,T y) { x = max(x,y); }
template <typename T> inline void chkmin(T &x,T y) { x = min(x,y); }
template <typename T> inline void read(T &x)
{
    T f = 1; x = 0;
    char c = getchar();
    for (; !isdigit(c); c = getchar()) if (c == '-') f = -f;
    for (; isdigit(c); c = getchar()) x = (x << 3) + (x << 1) + c - '0';
    x *= f;
}
inline int exp_mod(int a , int n)
{
        int b = a , res = 1;
        while (n > 0)
        {
                if (n & 1) res = 1LL * res * b % P;
                b = 1LL * b * b % P;
                n >>= 1;
        }
        return res;
}
inline int calc(int N , int M)
{
        vector< int > pr;
        int tmp = N / M;
        while (tmp != 1)
        {
                pr.push_back(f[tmp]);
                tmp /= f[tmp];        
        }    
        int sz = unique(pr.begin() , pr.end()) - pr.begin(); 
        int limit = m / M , res = 0;
        for (int i = 0; i < (1 << sz); ++i)
        {
                int sign = 1 , val = 1;
                for (int j = 0; j < sz; ++j)
                {
                        if (i & (1 << j))
                        {
                                sign *= -1;
                                val *= pr[j];
                        } 
                }
                res += 1LL * sign * (limit / val);
        }
        return res;
}

int main()
{
        
        read(m);
        for (rint i = 1; i <= m; ++i)
        {
                for (rint j = i; j <= m; j += i)
                {
                        d[j].push_back(i);                
                }
        }
        for (rint i = 2; i <= m; ++i)
        {
                if (!f[i])
                {
                        prime[++tot] = i;
                        f[i] = i;
                }
                for (int j = 1; j <= tot; ++j)
                {
                        int tmp = i * prime[j];
                        if (tmp >= MAXP) break;
                        f[tmp] = prime[j];
                        if (prime[j] == f[i]) break;
                } 
        }
        dp[1] = 1;
        for (rint i = 2; i <= m; ++i)
        {
                int res = 0;
                for (unsigned j = 0; j < d[i].size(); ++j)
                {
                        int D = d[i][j];
                        if (D != i) res = (res + 1LL * calc(i , D) * dp[D] % P) % P;
                }
                res = 1LL * res * exp_mod(m , P - 2) % P;
                res = (res + 1) % P;
                int fm = m - calc(i , i);
                res = 1LL * res * exp_mod(fm , P - 2) % P * m % P;
                dp[i] = res;
        }
        int ans = 0;
        for (rint i = 1; i <= m; ++i) ans = (ans + 1LL * exp_mod(m , P - 2) * dp[i] % P) % P;
        printf("%d
" , ans);
        
        return 0;
    
}