牛客练习赛84 F 牛客推荐系统开发之下班

好像只有我的公式是长这样的，但是跑的最快。
考虑到
其中，*是狄利克雷卷积，I是单位1。
进行杜教筛即可。
这里提一点，一个函数和mu卷是可以O(n)的求的。

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <map>
using namespace std;

typedef long long ll;

const int N = 5e6 + 100;
const int M = N - 100;
const int MOD = 1e9 + 9;

int tot;
int pri[N], f[N];
bool is_prime[N];

void upd(int &a, int b) {
    a += b; if (a >= MOD) a -= MOD;
}

void get_phi() {
    memset(is_prime, true, sizeof(is_prime));
    for (int i = 2; i < N; i++) {
        if (is_prime[i]) pri[tot++] = i;
        for (int j = 0; j < tot && i * pri[j] < N; j++) {
            is_prime[i * pri[j]] = false;
            if (i % pri[j] == 0) break;
        }
    }
    f[1] = 1;
    for (int i = 2; i <= M; i++) {
        f[i] = f[i - 1] + f[i - 2]; 
        if (f[i] >= MOD) f[i] -= MOD;
    }
    for (int i = 0; i < tot; i++) {
        for (int j = M / pri[i]; j; j--) {
            upd(f[j * pri[i]], MOD - f[j]);
        }
    }
    for (int i = 1; i <= M; i++) upd(f[i], f[i - 1]);
}

map<int, int> sf;

int fu[5][5], co[5][5], res[5][5];

void cal(int a[][5], int b[][5]) {
    for (int i = 0; i < 2; i++)
        for (int j = 0; j < 2; j++)
            for (int k = 0; k < 2; k++)
                fu[i][j] = (fu[i][j] + 1LL * a[i][k] * b[k][j]) % MOD;
    for (int i = 0; i < 2; i++) {
        for (int j = 0; j < 2; j++) {
            a[i][j] = fu[i][j];
            fu[i][j] = 0;
        }
    }
}

ll get(ll n) {
    n--;
    res[0][1] = res[1][0] = co[1][1] = 0;
    res[0][0] = res[1][1] = co[0][0] = co[0][1] = co[1][0] = 1;
    while (n > 0) {
        if (n & 1) cal(res, co);
        cal(co, co); n /= 2;
    }
    return res[0][0];
}

ll F(ll n) {
    if (n <= M) return f[n];
    if (sf.count(n)) return sf[n];
    int ans = get(n + 2) - 1;
    for (int l = 2, r; l <= n; l = r + 1) {
        r = n / (n / l);
        ans = (ans - (r - l + 1LL) * F(n / l)) % MOD;
    }
    ans = (ans % MOD + MOD) % MOD;
    return sf[n] = ans;
}

ll qpow(ll x, ll n) {
    ll res = 1;
    while (n > 0) {
        if (n & 1) res = res * x % MOD;
        x = x * x % MOD;
        n /= 2;
    }
    return res;
}

int n, k;
int a, b = 1;

int main() {
    //freopen("0.txt", "r", stdin);
    get_phi(); 
    scanf("%d%d", &n, &k);
    ll ans = 0;
    for (int l = 1, r; l <= n; l = r + 1) {
        r = n / (n / l);
        ll co = qpow(n / l, k);
        ans = (ans + co * (F(r) - F(l - 1))) % MOD;
    }
    ans = (ans % MOD + MOD) % MOD;
    printf("%lld\n", ans);
    return 0;
}