万能欧几里得入门学习笔记

• 前情提要：本人没有学习类欧几里得算法，大家不要学习我，该学的还是要学的。

给出一条直线 \(y = \dfrac {px + r} q\)，给出一个正整数 \(n\)，然后依次考虑每条横线 \(y = 1, 2, \dots, \lfloor \dfrac {pn + r} q \rfloor\) 和纵线 \(x = 1, 2, \dots, n\)。初始化一个矩阵为 \(I\)，从这条直线开始出发，若碰到一条横线则乘上矩阵 \(U\)，若碰到一条纵线则令乘上矩阵 \(R\)，若此时在一个整点上（即同时碰到一条横线和一条纵线），则先乘上 \(U\) 再乘上 \(R\)。

其中 \(R, U\) 为常量，我们的任务即是求出最终的矩阵。可以设计一个函数 \(\text{solve} (p, q, r, n, R, U)\)，返回值为最终的矩阵。
我们不用几何来表述该问题：有 \(n\) 个 \(R\)，第 \(i\) 个 \(R\) 之前需要出现 \(\lfloor \dfrac {pi + r} q \rfloor\) 个 \(U\)，且最后一个矩阵是 \(R\)，求所有矩阵的乘积。

我们可以尝试不断缩小问题规模：

• 如果 \(r\ge q\)，则开头一定会先出现 \(\lfloor \dfrac r q \rfloor\) 个 U，可以先去掉进行递归：\(\text{solve} (p, q, r, n, R, U) = U ^ {\lfloor \frac rq \rfloor} * \text{solve} (p, q, r \bmod q, n, R, U)\)。

• 如果 \(p\ge q\)，则每个 \(R\) 之前一定会先出现 \(\lfloor \dfrac p q \rfloor\) 个 U，同样可以递归：\(\text{solve} (p, q, r, n, R, U) = \text{solve} (p \bmod q, q, r, n, U ^ {\lfloor \frac p q \rfloor} R, U)\)。

• 此时 \(R\) 的数量已经比 \(U\) 的数量多了，可以考虑把整个序列进行 flip 进行递归。

设 \(m = \lfloor \dfrac {pn + r} q \rfloor\) 即 \(U\) 的个数。

如果 \(m = 0\) 那么直接返回 \(R ^ n\) 即可。

由 \(y = \dfrac {px + r} q \Rightarrow x = \dfrac {yq - r} p\)，第 \(i\) 个 \(U\) 后面第一个 \(R\) 的编号为 \(\lceil \dfrac {iq - r} p \rceil\)，则其前面有 \(\lfloor \dfrac {iq - r - 1} p \rfloor\) 个 \(R\)，那么感受一下认为我们需要递归 \(\text{solve} (q, p, -r - 1, m, U, R)\)。

但其实还有一些问题：

1. \(-r - 1\) 是负数，不能接受。

2. 最后一个矩阵不是 \(U\)。

对于第一个问题，我们考虑提出第一个 \(U\) 以及前面的 \(\lfloor \dfrac {q - r - 1} p \rfloor\) 个 \(R\) 并从原序列中删除，此时原来第 \(x + 1\) 个 \(U\) 现在变成了第 \(x\) 个 \(U\)，其前面 \(R\) 的个数为：\(\lfloor \dfrac {q(x + 1) - r - 1} p \rfloor - \lfloor \dfrac {q - r - 1} p \rfloor = \lfloor \dfrac {qx + (q - r - 1) \bmod p} p \rfloor\)，此时递归后变成了 \(r' = (q - r - 1) \bmod p\)，可以接受。

对于第二个问题，我们容易想到提出最后一个 \(U\) 后面的所有 \(R\) 并从原序列中删除，一共有 \(n - \lfloor \dfrac {qm - r - 1} p \rfloor\) 个 \(R\)。

所以有 \(\text{solve} (p, q, r, n, R, U) = R ^ {\lfloor \frac {q - r - 1} p \rfloor} * U * \text{solve} (q, p, (q - r - 1) \bmod p, m - 1, U, R) * R ^ {n - \lfloor \frac {qm - r - 1} p \rfloor}\)

同时也可以从几何的角度考虑这个问题，就是把原来这条直线进行翻折，然后删去左下角 \(x = 1\) 左边的部分并补齐，删去右上角多余的部分，再递归。

然后矩阵是能表示很多信息的，所以类欧能做的基本上万欧都能做，这两个东西似乎本质上是相同的。

类欧模板

点击查看代码

#include <bits/stdc++.h>
#define ll int
#define LL long long
#define uLL unsigned LL
#define fi first
#define se second
#define mkp make_pair
#define pir pair<ll, ll>
#define pb push_back
#define i128 __int128
using namespace std;
char buf[1 << 22], *p1, *p2;
// #define getchar() (p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, (1 << 22) - 10, stdin), p1 == p2)? EOF :
// *p1++)
template <class T>
const inline void rd(T &x) {
    char ch;
    bool neg = 0;
    while (!isdigit(ch = getchar()))
        if (ch == '-')
            neg = 1;
    x = ch - '0';
    while (isdigit(ch = getchar())) x = (x << 1) + (x << 3) + ch - '0';
    if (neg)
        x = -x;
}
const ll maxn = 310, mod = 998244353, iv = mod - mod / 2; const LL inf = 1e18;
ll power(ll a, ll b = mod - 2, ll p = mod) {
    ll s = 1;
    while (b) {
        if (b & 1)
            s = 1ll * s * a % p;
        a = 1ll * a * a % p, b >>= 1;
    }
    return s;
}
template <class T, class _T>
const inline ll pls(const T x, const _T y) { return x + y >= mod ? x + y - mod : x + y; }
template <class T, class _T>
const inline ll mus(const T x, const _T y) { return x < y ? x + mod - y : x - y; }
template <class T, class _T>
const inline void add(T &x, const _T y) { x = x + y >= mod ? x + y - mod : x + y; }
template <class T, class _T>
const inline void sub(T &x, const _T y) { x = x < y ? x + mod - y : x - y; }
template <class T, class _T>
const inline void chkmax(T &x, const _T y) { x = x < y ? y : x; }
template <class T, class _T>
const inline void chkmin(T &x, const _T y) { x = x < y ? x : y; }

struct Data {
    ll c, f, s1, s2, s3;
    const Data operator * (const Data t) const {
        return (Data) {pls(c, t.c), pls(f, t.f), (s1 + t.s1 + 1ll * c * t.f) %mod,
         (s2 + t.s2 + 1ll * c * c %mod * t.f + 2ll * c * t.s1) %mod,
         (s3 + t.s3 + (2ll * f + t.f + 1) * iv %mod * t.f %mod * c + 1ll * f * t.s1) %mod};
    }
};

Data qpow(Data a, ll b) {
    if(!b) return {};
    Data c = a; --b;
    while(b) {
        if(b & 1) c = c * a;
        a = a * a, b >>= 1;
    } return c;
}

Data solve(ll p, ll q, ll r, ll n, Data R, Data U) {
    if(!n) return {};
    if(r >= q) return qpow(U, r / q) * solve(p, q, r % q, n, R, U);
    if(p >= q) return solve(p % q, q, r, n, qpow(U, p / q) * R, U);
    ll m = (1ll * p * n + r) / q;
    if(!m) return qpow(R, n);
    return qpow(R, (q - r - 1) / p) * U * solve(q, p, (q - r - 1) % p, m - 1, U, R)
     * qpow(R, n - (1ll * q * m - r - 1) / p);
}

ll T, n, a, b, c;

void solve() {
    rd(n), rd(a), rd(b), rd(c);
    Data res = solve(a, c, b, n, {0, 1, 0, 0, 0}, {1, 0, 0, 0, 0});
    ll res1 = (res.s1 + b / c) %mod, res2 = (res.s2
     + 1ll * b / c * (b / c)) %mod, res3 = res.s3;
    printf("%d %d %d\n", res1, res2, res3);
}

int main() {
    rd(T); while(T--) solve();
    return 0;
}