BZOJ 4833: [Lydsy1704月赛]最小公倍佩尔数（数论 + 最值反演）

题面

令 ({(1+sqrt 2)}^n=e(n)+f(n)*sqrt2) ，其中 (e(n),f(n)) 都是整数，显然有 ({(1-sqrt 2)}^n=e(n)-f(n)*sqrt 2) 。

令 (g(n)) 表示 (f(1),f(2)…f(n)) 的最小公倍数，给定两个正整数 (n) 和 (p) ,其中 (p) 是质数，并且保证 (f(1),f(2)…f(n)) 在模 (p) 意义下均不为 (0) .

请计算(displaystyle sum_{i=1}^{n} i*g(i)) .其在模 (p) 的值。

(T le 210, sum n le 3 imes 10^6)

题解

其实它就是给你了 (f(n)) 的特征方程，其实可以解出来，具体来说 递推式到通相式 直接待定系数，解回去就把系数代回去就行了。

但是没有必要，(f(n)) 是个能打表的东西。

[f(n) = egin{cases} 0, &n = 0 \ 1, &n = 1 \ 2f(n - 1) + f(n-2) &mathrm{otherwise} end{cases} ]

然后它具有和 (fib(i)) 一样的性质， (gcd (f(i), f(j)) = f(gcd(i, j))) 。

至于证明，可以类似 (fib(i)) 的性质证明。

然后考虑最值反演：

不会的话可以看我 这篇博客 。

然后我们可以把这个对于 (mathrm{lcm}) 放到指数上来就是所有指数的 (max) ，那么我们可以用最值反演变成 (min) ，也就是 (gcd) 。

形式如下

[mathrm{lcm}(S) = prod_{T subseteq S, T ot = varnothing} gcd(T)^{(-1)^{|T|-1}} ]

然后我们利用这个性质进行推导。

[egin{align} g(n) &= prod_{T subseteq S, T ot = varnothing} gcd_{x in T}{f(x)}^{(-1)^{|T|+1}} \ &= prod_{T subseteq S, T ot = varnothing} f(gcd(T))^{(-1)^{|T|+1}}\ end{align} ]

然后考虑构造一个函数 (h) 满足

[f(n) = prod_{d|n} h(d) ]

由于在模素数意义下，所以必定存在 (h) 满足要求。

然后我们考虑每个 (h(d)) 的贡献：

[egin{align} Rightarrow g(n) &= prod_{d = 1}^{n} h(d)^{sum_{Tsubseteq S, T ot = varnothing}[d | gcd(T)](-1)^{|T|+1}} \ &= prod_{d = 1}^{n} h(d) end{align} ]

至于上面那个系数为什么 (=1) 呢？可以考虑 (T) 集合中所有为 (d) 倍数的数构成的所有集合，刚好是 (2^k - 1) 的，除去空集，那么剩下的一一抵消，剩下一个 (1) 的贡献。

考虑求出 (h) 是可以枚举倍数来优化到 (O(n ln n)) 的。

具体来说假设我们 (<n) 的所有 (h) 都求好了，那么当前的 (displaystyle h(n) = f(n) prod_{d|n, d ot = n} frac{1}{h(d)}) 。

其实可以用莫比乌斯反演推回去的，但没必要筛 (mu) 。

那么能不能用线性筛优化到 (O(n)) 呢？似乎是不行的，因为 (h) 不是积性的，我们知道 (h(i)) 与 (h(p)) 无法推出 (h(i imes p)) 。

然后最后复杂度就是 (O(sum n ln n)) 的。

总结

对于有些看起来不太好求的 (mathrm{lcm}) 题，考虑最值反演。然后考虑构造积性函数，最后求出积性函数求回来就行了。

代码

/**************************************************************
    Problem: 4833
    User: zjp_shadow
    Language: C++
    Result: Accepted
    Time:5908 ms
    Memory:13024 kb
****************************************************************/
 
#include <bits/stdc++.h>
 
#define For(i, l, r) for(register int i = (l), i##end = (int)(r); i <= i##end; ++i)
#define Fordown(i, r, l) for(register int i = (r), i##end = (int)(l); i >= i##end; --i)
#define Set(a, v) memset(a, v, sizeof(a))
#define Cpy(a, b) memcpy(a, b, sizeof(a))
#define debug(x) cout << #x << ": " << (x) << endl
#define DEBUG(...) fprintf(stderr, __VA_ARGS__)
 
using namespace std;
 
typedef long long ll;
 
template<typename T> inline bool chkmin(T &a, T b) {return b < a ? a = b, 1 : 0;}
template<typename T> inline bool chkmax(T &a, T b) {return b > a ? a = b, 1 : 0;}
 
inline int read() {
    int x(0), sgn(1); char ch(getchar());
    for (; !isdigit(ch); ch = getchar()) if (ch == '-') sgn = -1;
    for (; isdigit(ch); ch = getchar()) x = (x * 10) + (ch ^ 48);
    return x * sgn;
}
 
void File() {
#ifdef zjp_shadow
    freopen ("4883.in", "r", stdin);
    freopen ("4883.out", "w", stdout);
#endif
}
 
const int N = 1e6 + 1e3;
 
int Mod, h[N], g[N], ans[N];
 
ll fpm(ll x, int power) {
    ll res = 1;
    for (; power; power >>= 1, (x *= x) %= Mod)
        if (power & 1) (res *= x) %= Mod;
    return res;
}
 
void Solve(int maxn) {
    h[1] = 1;
    For (i, 2, maxn)
        h[i] = (h[i - 1] * 2ll + h[i - 2]) % Mod;
    For (i, 2, maxn) {
        int inv = fpm(h[i], Mod - 2);
        for (int j = i * 2; j <= maxn; j += i)
            h[j] = 1ll * h[j] * inv % Mod;
    }
    g[0] = 1; int ans = 0;
    For (i, 1, maxn) {
        g[i] = 1ll * g[i - 1] * h[i] % Mod;
        ans = (ans + 1ll * g[i] * i) % Mod;
    }
    printf ("%d
", ans);
}
 
int main () {
     
    File();
 
    int cases = read();
    while (cases --) {
        int n = read(); Mod = read(); Solve(n);
    }
 
    return 0;
 
}