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 \times 10^6\)
题解
其实它就是给你了 \(f(n)\) 的特征方程，其实可以解出来，具体来说 递推式到通相式 直接待定系数，解回去就把系数代回去就行了。
但是没有必要，\(f(n)\) 是个能打表的东西。
然后它具有和 \(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 \not = \varnothing} \gcd(T)^{(-1)^{|T|-1}} \]
然后我们利用这个性质进行推导。
然后考虑构造一个函数 \(h\) 满足
由于在模素数意义下，所以必定存在 \(h\) 满足要求。
然后我们考虑每个 \(h(d)\) 的贡献：
至于上面那个系数为什么 \(=1\) 呢？可以考虑 \(T\) 集合中所有为 \(d\) 倍数的数构成的所有集合，刚好是 \(2^k - 1\) 的，除去空集，那么剩下的一一抵消，剩下一个 \(1\) 的贡献。
考虑求出 \(h\) 是可以枚举倍数来优化到 \(O(n \ln n)\) 的。
具体来说假设我们 \(<n\) 的所有 \(h\) 都求好了，那么当前的 \(\displaystyle h(n) = f(n) \prod_{d|n, d \not = n} \frac{1}{h(d)}\) 。
其实可以用莫比乌斯反演推回去的，但没必要筛 \(\mu\) 。
那么能不能用线性筛优化到 \(O(n)\) 呢？似乎是不行的，因为 \(h\) 不是积性的，我们知道 \(h(i)\) 与 \(h(p)\) 无法推出 \(h(i \times 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\n", ans);
}
int main () {
File();
int cases = read();
while (cases --) {
int n = read(); Mod = read(); Solve(n);
}
return 0;
}
