// #pragma GCC optimize(3, "Ofast", "inline", "unroll-loops")
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 1100010;
const int mod = 998244353;
const int inf = 1e18;
using ld = long double;
using ull = unsigned long long;
using i128 = __int128;
const ull base = 13331;
namespace Luminescent
{
const double pi = acos(-1);
const ld pi_l = acosl(-1);
struct DSU
{
int fa[N];
inline DSU() { iota(fa, fa + N, 0); }
inline void init(int maxn) { iota(fa, fa + maxn + 1, 0); }
inline int find(int x) { return x == fa[x] ? x : fa[x] = find(fa[x]); }
inline int merge(int x, int y)
{
x = find(x), y = find(y);
if (x != y)
return fa[x] = y, 1;
return 0;
}
};
inline void add(int &x, int a) { x = (x + a) % mod; }
inline void sub(int &x, int a) { x = (x - a + mod) % mod; }
inline int power(int a, int b, int c)
{
int sum = 1;
while (b)
{
if (b & 1)
sum = 1ll * sum * a % c;
a = 1ll * a * a % c, b >>= 1;
}
return sum;
}
inline int inversion(int x) { return power(x, mod - 2, mod); }
inline int inversion(int x, int mod) { return power(x, mod - 2, mod); }
inline int varphi(int x)
{
int phi = 1;
for (int i = 2; i * i <= x; ++i)
if (x % i == 0)
{
phi *= (i - 1);
x /= i;
while (x % i == 0)
phi *= i, x /= i;
}
if (x > 1)
phi *= (x - 1);
return phi;
}
}
using namespace Luminescent;
namespace Poly
{
const int g = 3;
int rev[N];
void ntt(int *a, int n, int mode)
{
int bit = 1;
while ((1 << bit) < n)
++bit;
for (int i = 0; i < n; ++i)
{
rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (bit - 1));
if (i < rev[i])
swap(a[i], a[rev[i]]);
}
for (int l = 2; l <= n; l <<= 1)
{
int x = power(g, (mod - 1) / l, mod);
if (mode == 1)
x = inversion(x);
for (int i = 0; i < n; i += l)
{
int v = 1;
for (int j = 0; j < l / 2; ++j, v = v * x % mod)
{
int v1 = a[i + j], v2 = a[i + j + l / 2] * v % mod;
a[i + j] = (v1 + v2) % mod, a[i + j + l / 2] = (v1 - v2 + mod) % mod;
}
}
}
}
// calc convolution: c[i] = \sum\limits_{j=0}^i (a[j] * b[i - j])
void convolution(int *a, int n, int *b, int m, int *c)
{
int tn = n, tm = m;
n = n + m + 2;
while (__builtin_popcount(n) > 1)
++n;
// cerr << "n = " << n << '\n';
for (int i = tn + 1; i <= n + 1; ++i)
a[i] = 0;
for (int i = tm + 1; i <= n + 1; ++i)
b[i] = 0;
ntt(a, n, 0), ntt(b, n, 0);
for (int i = 0; i < n; ++i)
c[i] = a[i] * b[i] % mod;
ntt(c, n, 1);
const int inv_n = inversion(n);
for (int i = 0; i <= n + m; ++i)
c[i] = c[i] * inv_n % mod;
}
}
namespace Loyalty
{
inline void init() { }
int y[N], n, m;
int fac[N], inv[N], ifac[N];
int A[N], B[N], C[N], P[N];
inline void main([[maybe_unused]] int _ca, [[maybe_unused]] int atc)
{
cin >> n >> m;
for (int i = 0; i < 2; ++i)
fac[i] = inv[i] = ifac[i] = 1;
for (int i = 2; i < N; ++i)
{
fac[i] = fac[i - 1] * i % mod;
inv[i] = mod - inv[mod % i] * (mod / i) % mod;
ifac[i] = ifac[i - 1] * inv[i] % mod;
}
for (int i = 0; i <= n; ++i)
cin >> y[i];
auto coef = [&](int x) { return (x & 1) ? (mod - 1) : 1; };
for (int i = 0; i <= n; ++i)
A[i] = ifac[i] * y[i] % mod * coef(n - i) % mod * ifac[n - i] % mod;
for (int i = 0; i <= n + n; ++i)
B[i] = inversion(m + i - n);
Poly::convolution(A, n + n, B, n + n, C);
int fac_m = 1, ifac_m = 1;
for (int i = 2; i <= m; ++i)
fac_m = fac_m * i % mod;
for (int i = 2; i <= m - n - 1; ++i)
ifac_m = ifac_m * i % mod;
ifac_m = inversion(ifac_m);
for (int k = 0; k <= n; ++k)
{
P[k] = fac_m * ifac_m % mod;
fac_m = fac_m * (m + k + 1) % mod;
ifac_m = ifac_m * inversion(m + k - n) % mod;
}
for (int k = 0; k <= n; ++k)
cout << C[n + k] * P[k] % mod << ' ';
cout << '\n';
}
}
signed main()
{
// freopen("1.in", "r", stdin);
// freopen("1.out", "w", stdout);
cin.tie(0)->sync_with_stdio(false);
cout << fixed << setprecision(15);
int T = 1;
// cin >> T;
Loyalty::init();
for (int ca = 1; ca <= T; ++ca)
Loyalty::main(ca, T);
return 0;
}
:::
002. 代码源 NOIp 模拟赛 Day3 宇宙（universe）
这个题真难吧（）
首先可以观察到每一个时刻在同一个维度上加多次坐标是没啥用的，于是容易想到看增加的距离的和是否足够在当前时刻逃离黑洞。然后可以发现这个东西具有单调性，所以考虑先二分时间 tim 判断该时刻是否仍能不被黑洞吃掉。然后内层需要增加的距离的和形如 \sum\limits_{i=1}^n\max(0,tim+1-a_i) 的形式，因此想到给 a 数组从小到大排序，然后二分出 \max 函数的转折点，维护前缀和算答案，总时间复杂度为 O(n\log^2n) 据说卡常之后可以过，但是我 T 了
然后注意到随着 k 增大，需要考虑的维度数量是单调递增的。因此想到用双指针维护这个变化过程，可以优化到 O(n\log n) 解决。
:::success[Code]
int a[N], pre[N];
inline void main([[maybe_unused]] int _ca, [[maybe_unused]] int atc)
{
int cccccccc, n;
cin >> cccccccc >> n;
for (int i = 1; i <= n; ++i)
cin >> a[i];
sort(a + 1, a + n + 1);
for (int i = 1; i <= n; ++i)
pre[i] = pre[i - 1] + a[i] - 1;
for (int k = 1, i = 1; k < n; ++k)
{
while (i < n && (i <= k || a[i + 1] * i - pre[i] <= k * a[i + 1]))
++i;
cout << pre[i] / (i - k) << ' ';
// cout << i << ' ';
}
cout << '\n';
}
:::
003. 代码源 NOIp 模拟赛 Day3 跳跃（jump）
严格还原了考场上的思考过程（）
注意到 C 性质给的分很多，于是考虑 C 性质有什么用。注意到此时保证串内不存在 k 个连续的 0，然后可以想到每存在 k 个连续的 0，都可以假装这些 0 不存在，然后到这里的时候就加一次跳跃再加一次经过 0，对串 s 重标号一下就可以转化为 C 性质的情况即不存在连续的 k 个 0。
判断当前二分到的总长度 x 是否合法是简单的，直接差分约束做即可。但是问题是若当前不合法，则无法知道合法区间是在 x 的左侧还是在 x 的右侧，这里考虑使用一些巧妙的手法：
对于建边关系而言：
对于 1 操作：
若有 a<b，则其距离起点 1 的数组 x 的关系形如 x_b-x_a\ge c，转化成基本模型是 x_a-x_b\le -c，该边边长和二分到的 x 没有任何关系。
若有 a>b，则其 x 值的关系形如 x_b-x_a+x\ge c，转化成基本模型即 x_b-x_a\ge c-x，x_a-x_b\le x-c，该边的边长和 x 成一倍增加关系。
对于 2 操作：
若有 a<b，则其距离起点 1 的数组 x 的关系形如 x_b-x_a\le c，该边边长和二分到的 x 没有任何关系。
若有 a>b，则其距离起点 1 的数组 x 的关系形如 x_b-x_a+x\le c 即 x_b-x_a\le c-x，该边边长和 x 成负一倍增加关系。
然后还有相邻两个元素之间的关系，对于 i=1\ldots n-1，其形如 x_i+1\le x_{i+1}，容易发现该边边长和二分到的 x 没有任何关系。而对于 i=n 的特殊情况，此时有 x_n\le x_1+x-1 即 x_n-x_1\le x-1，此时该边边长和 x 成一倍增加关系。
于是考虑在差分约束建边的时候，额外增加一个信息表示其和二分到的值 x 之间的关系，该关系的值一定是 -1,0,1 中的一个。然后若总长度为 x 是不合法的，则考虑找出其所在的负环（这个可以在 spfa 的过程中通过记录路径得到），容易证明若其负环上所有倍率边的和为 0，则此时不管怎么调整 x 的值都必然会存在负环，即必然无解。否则，若该值 >0，则通过把 x 的值往大调整就可以使得该负环上边权值的和增大直到其权值之和为正。
通过这个性质可以得出下面的结论：对于一张给定的竞赛图而言，竞赛图的 SCC 数量为 \sum\limits_{i=1}^n{[\sum\limits_{j=1}^id_j'=\binom i2]}，其中 d' 数组是图上所有点的入度 \deg 升序排序得到的新数组。
直接暴力枚举竞赛图上每条边的朝向，时间复杂度为 O(2^{n^2})，期望得分 0 分。
Algorithm 2
暴力枚举朝向还是太？了，有没有更正常的做法？
考虑枚举划分的位置 i（即对给定的 i 上面艾弗森括号内取等），计数该条件成立时的概率。考虑这个条件成立的另一个定义，其实就是要求该竞赛图上的某 i 个点组成的集合 L 连向另外 n-i 个点组成的集合 R 的边都是由 L 集合单向连向 R 集合的。很显然这个条件成立的概率就是 2^{-i(n-i)}。而从 n 个点中选 i 个点划分出 L,R 两个集合的方案数则是 \binom ni。
于是考虑从度数方面入手，动态维护若加入 / 删除当前区间内的边，森林是否仍然为 / 不为毛毛虫森林。直接加入一个点可能存在后效性（加入一个点 x，此时 x 必然为叶子结点，则其会影响和她相邻的唯一结点 y 的信息，而 y 又会影响到和其相邻的另一个点 w 的信息），维护起来比较麻烦。然后我编了一个奇奇怪怪的做法，又维护了另一个数组信息，表示每个点相邻的点在当前导出子图中的点编号的异或，根据异或的一些性质，可以快速的找出后效性产生的结点 w，而这个数组也同样是容易维护的。
题目给出了明确的提示“把 x_i,y_i 看成一个整体组来考虑”，于是容易想到 Lagrange 插值插出 f 函数，根据经典结论可知插出的函数恰有 n 项，然后再编一个大于 \max y_i 的模数（例如 10^9+7）在模意义下插值。此时 f 函数每一项的系数都不超过 2^{30}，直接转成二进制扔到 s 串里传输即可。
而对于逆变换，直接套 Lagrange 插值板子即可。
029. P6790 [SNOI2020] 生成树
前置知识：广义串并联图
核心操作是：
删一度点：若图上一个点 u 的度数为 1，则直接删除该点。
缩二度点：若图上一个点 u 只和 v,w 两个点相连，则删除 u 点，以及 u,v 边 u,w 边，连接 v,w 边。
一些比较难受的地方是可持久化线段树空间要开足（虽然这个题空间限制是 1G 但是还是很卡空间，不能全开 long long）。
2025.11.14
037. CF2085F2 Serval and Colorful Array (Hard Version)
先考虑比较暴力的做法。枚举选择的子序列最中间的位置 i，这样一来根据经典贪心结论，左右两侧都分别需要选 \frac k2 个（对 2\mid k 的情况需要处理左边多选一个还是右边多选一个）。对每个值不和 i 相等的值，求出 i 左右两侧第一个和其值相同的位置 l_j,r_j，贪心排序选取就可以做到 O(n^2\log n)。
这个元素和前面某个和其原来不在同一序列的元素匹配，此时可以枚举上一个匹配位置 k，转移方程形如 f_{i,j}\leftarrow f_{i-1,j-1}\times c_i。这个时候注意到转移式和 k 没有关系，因此只需要计数 i 前面和 i 原来不在同一序列的元素的数量即可，这个可以一个前缀和求出来。
最后将答案乘以 (n!)^{-1} 即可。总时间复杂度为 O(n^2)，可以通过该题。
042. 代码源 NOIp 模拟赛 Day6 诅咒
弱化版是 QOJ8723 2024北京市赛 乘二。
考虑优化。注意到该题中堆的操作都比较特殊，因此使用 [P2827 [NOIP 2016 提高组] 蚯蚓](https://www.luogu.com.cn/problem/P2827) 一题的套路，用两个有序队列 $q_1,q_2$ 来模拟出一个堆。具体细节见代码。
一些特殊情况是要特判 $k=1$ 和 $\min a=0$ 的 corner case。
:::success[AC 代码]
```cpp line-numbers
namespace Loyalty
{
int a[N];
queue<int> q1, q2;
inline void push(int x) { q2.emplace(x); }
inline int top()
{
int res = inf;
if (q1.size())
res = q1.front();
if (q2.size())
res = min(res, q2.front());
return res;
}
inline void pop()
{
int o = top();
if (q1.size() && q1.front() == o)
q1.pop();
else
q2.pop();
}
inline void init() {}
inline void main([[maybe_unused]] int _ca, [[maybe_unused]] int atc)
{
int n, m, k, s = 0;
cin >> n >> m >> k;
for (int i = 1; i <= n; ++i)
cin >> a[i];
if (!m)
{
int s = 0;
for (int i = 1; i <= n; ++i)
s = (s + a[i]) % mod;
cout << s << '\n';
return;
}
if (k == 1 || !a[1])
{
int res = 0;
sort(a + 1, a + n + 1);
for (int i = 1; i <= n; ++i)
res = (res + a[i]) % mod;
cout << res << '\n';
return;
}
sort(a + 1, a + n + 1);
for (int i = 1; i <= n; ++i)
q1.emplace(a[i]);
while (1)
{
int t = top();
if (t * k <= a[n])
{
pop();
push(t * k);
if (!--m)
break;
}
else
break;
}
// for (int i = 1; i <= n; ++i)
// a[i] = top(), pop();
// sort(a + 1, a + n + 1);
// for (int i = 1; i <= n; ++i)
// q1.emplace(a[i]);
for (int i = 1; i <= m % n; ++i)
{
int t = top();
pop();
push(t * k);
}
for (int i = 1; i <= n; ++i)
a[i] = top(), pop();
for (int i = 1; i <= n; ++i)
a[i] = a[i] % mod * power(k, m / n, mod) % mod, s = (s + a[i]) % mod;
cout << s << '\n';
}
}
```
:::
### 044. [代码源 NOIp 模拟赛 Day6 区间](http://oj.daimayuan.top/contest/426/problem/3444)
原题是 [P12448 [COTS 2025] 观草 / Trava](https://www.luogu.com.cn/problem/P12448)。
问题询问的式子是 $\sum\limits_{i=1}^{n-k+1}\max\limits_{i\le j\le i+k-1}a_j$，看上去十分不能维护。于是考虑对这个式子做一点小小的转化：
$$
\begin{aligned}
&\sum\limits_{i=1}^{n-k+1}\max\limits_{i\le j\le i+k-1}a_j\\
=&\sum\limits_{1\le p\le V}\sum\limits_{i=1}^{n-k+1}[\max\limits_{i\le j\le i+k-1}a_j\ge p]\\
=&(n-k+1)V-\sum\limits_{1\le p\le V}\sum\limits_{i=1}^{n-k+1}[\min\limits_{i\le j\le i+k-1}a_j<p]\\
\end{aligned}
$$
上面的分析套路的把一个序列转化成了 $01$ 序列。现在考虑在这里搞事情。设另一个序列 $b$ 满足 $b_i=[a_i<p]$（对每个 $p$ 而言 $b$ 序列都是相同的），此时上面艾弗森括号内的条件就变为了 $[\min\limits_{i\le j\le i+k-1}b_j]$ 即 $\forall j\in[i,i+k-1],b_j=\text{True}$。因此想到只考虑 $b$ 序列中值为 $1$ 的位置。设开头结尾和每一个 $b$ 序列中 $1$ 的位置将序列分隔为了 $x$ 段，长度分别为 $len_1,len_2,\ldots,len_x$，则答案即为：$(n-k+1)V-\sum\limits_{1\le p\le V}\sum\limits_{i=1}^x\max(0,len_i-k+1)$。
此时这个东西仍然是没法求的。考虑到 $len$ 数组的定义就是序列中的极长全 $0$ 段。考虑按照从小到大的顺序扫描 $p$，经典结论（也可以均摊分析）可以发现每个位置的值最多会变化一次，也就是说不同极长全 $0$ 段是 $O(n)$ 级别的，于是想到拿单调栈对每个 $i$ 维护出其为最值时的极长区间 $[l_i,r_i]$，然后做一个后缀修改即可。
然后你还需要维护一个形如 $\max(0,len_i-k+1)$ 的式子，这是 [P3586 [POI 2015 R2] 物流 Logistics](https://www.luogu.com.cn/problem/P3586) 的经典套路，可以维护两个树状数组然后做个差求和做到小常数单 $\log$ 求解。
单点 $+1$ 也是简单的，考虑到这样只会更改 $O(1)$ 个极长全 $0$ 段，因此想到直接左右两个线段树二分求出该元素对应的极长 $0$ 区间然后将其分割为两个极长全 $0$ 区间，分别在树状数组和线段树上做修改即可。总时间复杂度为 $O(n\log n)$，可以通过该题。
## 2025.11.16
### 045. [P11536 [NOISG 2023 Finals] Curtains](https://www.luogu.com.cn/problem/P11536)
根据 CF2163D1 的结论可知，对于一个左端点，只有右端点最大的区间是有用的，因此只需要保留最多 $n$ 个区间。
考虑对询问区间和给定的区间分别排序然后双指针扫描，对于一个询问区间，将右端点在其右端点左侧的覆盖区间加入到线段树中，则此时询问区间 $[ql,qr]$ 合法的充要条件是：
+ $[ql,qr]$ 中每个位置都被至少覆盖一层。
+ $[ql,qr]$ 中所有位置对应的覆盖左端点的最大值的 $\min$ 值 $=ql$。
发现这个东西就是区间取 $\max$ 区间查询 $\min$，可以简单线段树维护，总时间复杂度为 $O(n\log n)$ 可以通过该题。
### 046. [P10045 [CCPC 2023 北京市赛] 线段树](https://www.luogu.com.cn/problem/P10045)
没见过就不会做系列
因为 $a_i$ 时刻为奇数，所以考虑把 $a_i$ 表示成 $2b_i+1$ 的形式。此时区间 $[l,r]$ 的乘积即可以表示为 $\prod\limits_{i=l}^r(2b_i+1)$ 的形式。又因为答案是对 $2^{20}$ 取模的，且有 $2\mid 2b_i$，因此上面的这个式子只用考虑有最多 $19$ 个选 $2b_i$ 的情况。
考虑直接用线段树维护这个东西。记 $cnt_i$ 表示在当前线段树结点对应区间 $[l,r]$ 上选 $i$ 个 $2b_i$ 作为乘积的方案数。合并枚举左右两侧选择的 $2b_i$ 数量然后简单计算即可。根据上面的分析只需要维护 $i=0\ldots 19$ 时的 $cnt_i$ 值即可。
总时间复杂度为 $O(nk^2\log n)$，其中 $k$ 取 $20$，卡常后可以过。
一个卡常技巧是答案用 `unsigned int` 存，众所周知 `unsigned int` 的乘法取模比 `int` 快，而且这个东西对 $2^{32}$ 自然溢出，然后就可以卡着时间限制通过。谴责出题人卡常。
### 047. [[CXOI R5-A] Decomposition](https://www.luogu.com.cn/problem/T676772?contestId=290075)
计算初始答案是简单的，发现题目问的东西就是对每个结点 $i$，求 $1\sim i$ 路径上的左链数量。考虑一次树上 dp 求出 $siz_i$ 表示 $i$ 点为根的子树的结点数量，$f_i$ 表示 $1\sim i$ 路径上的左链数量，这两个信息都是可以 $O(n)$ 求出的。因此初始答案即为 $\sum\limits_{i=1}^nf_i$。
然后还有 $k$ 次查询操作，考虑经典 exchange-argument 技巧，对一个同时拥有左右儿子的结点 $u$，考虑计算若交换两个儿子，则对答案会额外产生多少贡献（$\Delta$）。设其左右两个儿子的子树大小分别为 $s_l,s_r$，则有 $\Delta=s_l-s_r$。
把所有 $\Delta$ 排序然后贪心选前 $k$ 大且 $\Delta$ 为正数的结点交换左右儿子即可，总时间复杂度为 $O(n\log n)$ 可以通过该题。
### 048. [[CXOI R5-B] Calculator](https://www.luogu.com.cn/problem/T683746?contestId=290075)
你一定会用牛顿迭代法开根。
### 049. [P14518 [NFLSPC #8] APLSPC](https://www.luogu.com.cn/problem/P14518)
怎么还有 Quine，并非 OI 题，就不写题解了，可以去看第一篇题解。
## 2025.11.17
### 050. [P7962 [NOIP2021] 方差](https://www.luogu.com.cn/problem/P7962)
考虑沿用 CBC022E 的套路，把题面中给的操作扔到 $a$ 的差分数组里，容易发现其实就是交换两个相邻的数的差分值。
套路的把方差乘以 $n^2$ 拆成 $n\sum\limits_{i=1}^na_i^2-(\sum\limits_{i=1}^na_i)^2$ 的形式，然后考虑 dp 处理。设 $f_{i,j}$ 表示当前处理完了前 $i-1$ 个差分值，此时平方项的和值为 $j$，最少花费是多少。把方差的形式用差分的前缀和的形式展开然后合并同类项，可以发现最后得到的式子是一个单峰的形式，因此对于前面给出的 dp，插入一个新的元素要么插入到其最左侧要么插入到其最右侧。转移是容易的。
直接这样做时间复杂度为 $O(n^2a)$，难以通过。考虑继续优化。注意到 dp 上值不为 $0$ 的数是很少的，所以考虑只转移 dp 不为 $0$ 的情况，容易计算得出这样的位置的数量是只有 $O(a)$ 个的。因此此时总时间复杂度为 $O(na^2)$，可以通过该题。
### 051. [P3430 [POI 2005] DWU-Double-row](https://www.luogu.com.cn/problem/P3430)
感觉是挺套路的题，考虑相邻位置之间建边，若值为 $i$ 的两个位置分别出现在第 $x$ 列和第 $y$ 列，那么：
+ 若两个位置出现在同一行，则必须换，在 $x,y$ 之间连一条权值为 $1$ 的双向边。
+ 若两个位置不出现在同一行，那么必须相对不换，在 $x,y$ 之间连一条权值为 $0$ 的双向边。
然后容易发现这是个二分图。直接二分图染色，对每个连通块把出现次数少的那个颜色对应的所有位置全部换位置即可满足条件。
总时间复杂度为 $O(n)$，可以通过该题。
### 052. [QOJ 5 在线 $O(1)$ 逆元](https://qoj.ac/problem/5)
科技的力量.jpg
前置知识：[Farey 序列](https://oi-wiki.org/math/number-theory/stern-brocot/)以及其相关理论，根号平衡。
这里稍微提一下什么是 Farey 序列 以及其性质：
+ $F_n$ 是第 $n$ 阶的 Farey 序列，序列中存储所有不可约的分母 $\le n$ 的值在 $[0,1]$ 之间的分数（这里认为 $\frac01,\frac11$ 也是不可约的分数）按照其值从小到大排序得到的有序分数序列。
+ 性质 $1$：对于 $F_{m-1}$ 序列上任意两个相邻分数 $\frac ab,\frac cd$，在 $F_m$ 序列上一定**按照顺序连续的**存在分数 $\frac ab,\frac{a+c}{b+d},\frac bd$。
+ 推论 $1$：对于 $F_m$ 中任意两个相邻分数 $\frac ab,\frac cd$，必然有 $bc-ad=1$。
+ 性质 $2$：$|F_n|=1+\sum\limits_{i=1}^n\varphi(i)$。
+ 性质 $3$：对任意整数 $n\ge 2$ 和任意实数 $v\in[0,1]$，一定可以在 $F_{n-1}$ 序列中找到至少一个分数 $\frac xy$，满足 $|v-\frac xy|\le \frac1{yn}$。
+ 推论 $3$：分数 $\frac xy$ 一定是数 $v$ 在序列 $F_{n-1}$ 中向左查或者向右查能找到的第一个分数。
+ 性质 $4$：对于 $F_n$ 序列内的每个分数 $\frac xy$，$\lfloor\frac{xn^2}{2y}\rfloor$ 的值互不相同。
回到这个题目。考虑固定 $n$，记 $v=\frac ap$。根据上面得到的性质和推论，可以知道此时必然存在一个分数 $\frac xy\in F_{n-1}$ 满足 $|\frac ap-\frac xy|\le\frac1{yn}$。此时把不等式的左右两侧同时乘以 $py$（显然有 $py>0$），得到：$|ay-px|\le\lfloor\frac pn\rfloor$。
记 $u=ay-px$，则此时又能得出两个结论：
+ $ay\equiv u\pmod p
考虑把上面提到的结论 1 改写，得到：y\equiv u\times a^{-1}\pmod p。因为这里想要求的是 a 在模 p 意义下的逆元即 a^{-1}，因此想到在同余式两侧同时乘以 u^{-1}，得到：a^{-1}\equiv u^{-1}\times y\pmod p。
然后再利用结论 2 即 |u|\le\lfloor\frac pn\rfloor 也就是 u 的绝对值很小，因此对这样的 a 只需在对应 Farey 序列上找出其对应的分数 \frac xy，计算 u=ay-px，即可套用公式 a^{-1}\equiv u^{-1}y\pmod p 解决。而此时 u^{-1} 可以直接线性求逆元。
然而此时存在 u<0 的情况。对这个东西的处理是比较容易的，有 p^{-1}\equiv -(-p)^{-1}\pmod p 然后转化成上一种情况了。
#include "inv.h"
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 2100010;
int pre[N], idx, mod, n, m, Inv[N];
pair<int, int> frac[N], farey[N];
#undef int
void init(int p)
#define int long long
{
n = cbrt(p), m = n * n;
for (int i = 1; i < n; ++i)
for (int j = 0; j <= i; ++j)
{
int val = j * m / i;
if (!pre[val])
pre[val] = 1, frac[val] = {j, i};
}
for (int i = 1; i <= m; ++i)
pre[i] += pre[i - 1];
idx = 0;
for (int i = 0; i <= m; ++i)
if (frac[i].first || frac[i].second)
farey[idx++] = frac[i];
Inv[0] = Inv[1] = 1;
for (int i = 2; i <= p / n; ++i)
Inv[i] = Inv[p % i] * (p - p / i) % p;
mod = p;
}
#undef int
int inv(int a)
#define int long long
{
int val = a * m / mod;
if (pre[val] < idx)
{
auto [x, y] = farey[pre[val]];
int w = a * y - x * mod;
if (abs(w) <= mod / n)
{
int inv_t;
if (w >= 0)
inv_t = Inv[w];
else
inv_t = (mod - Inv[-w]) % mod;
return y * inv_t % mod;
}
}
if (pre[val] - 1 >= 0)
{
auto [x, y] = farey[pre[val] - 1];
int w = a * y - x * mod;
if (abs(w) <= mod / n)
{
int inv_t;
if (w >= 0)
inv_t = Inv[w];
else
inv_t = (mod - Inv[-w]) % mod;
return y * inv_t % mod;
}
}
throw;
return -1;
}
考虑证明该结论。上面这个东西是 \mu 函数和 d 函数的狄利克雷卷积，众所周知 \mu 和 d 都是积性函数，所以 f(n)=\mu * d 也同样是积性函数。根据经典套路可知，积性函数只需要求出 f(p)（p 是质数）的值就可以直接推出其他所有值（例题：CBC013E，给定 k 求所有 i^k 的值，f(x)=x^k 是完全积性函数，可以通过筛出质数位置的值然后递推出来）
考虑从 g 入手推出 f 的式子，此时有 g(p^k)=(p^k)^2=p^{2k}=\sum\limits_{i=0}^kf(p^i)，移项可得 f(p^i)=p^{2k}-\sum\limits_{i=0}^{k-1}f(p^i)=p^{2k}-g(p^{k-1})=p^{2k}-p^{2k-2}。