B3646

· · 题解

数列前缀和 3

这是一道洛谷夏令营的作业题,出这道题的目的是强调非交换群的区间和必须注意左右乘的区别。

考虑维护矩阵的前缀积:s_i 表示前 i 个矩阵的乘积。 注意对于一个 [l, r] 的询问,\prod\limits_{i = l}^r a_i \neq s_r \times s_{l-1}^{-1} 。例如:CD \neq ABCD \times (AB)^{-1}

但是注意到 (AB)^{-1} = B^{-1} A^{-1}。我们只要将 (AB)^{-1} 左乘上 ABCD,就有 (AB)^{-1} \times ABCD = (B^{-1}(A^{-1}A)B)CD = CD

所以我们有 \prod\limits_{i = l}^r a_i = s_{l - 1}^{-1} \times s_r。维护前缀积,做一个矩阵求逆即可。时间复杂度 O((n+q) k^3)

矩阵求逆可以用高斯消元法,但因为 k 只有 23,所以也可以手搓一下伴随矩阵,下面介绍一下伴随矩阵法求逆矩阵。

n 阶矩阵 A=(a_{ij})_{n \times n},划去第 i 行第 j 列后剩余的 (n-1) 阶行列式的值称为元素 a_{i,j} 的余子式,记为 M_{i, j}。记 A_{i, j}=(-1)^{i+j}M_{i,j}a_{i, j} 的代数余子式。

对于二阶矩阵 $A = \begin{pmatrix} a & b \\ c & d\end{pmatrix}$,直接套用上述结论得到 $A^{-1} = \frac{1}{ad - bc}\begin{pmatrix} d & -b \\ -c & a\end{pmatrix}$;对于三阶矩阵,其行列式只有六项,所有的代数余子式均为二阶行列式,都可以轻松算出。 ```cpp #include <array> #include <iostream> #include <algorithm> const int p = 1145141; typedef long long int ll; int n, k, q; std::array<int, p> inv; struct Matrix { ll A[3][3]; inline Matrix operator*(const Matrix& o) const { Matrix ret; for (int i = 0; i < k; ++i) { for (int j = 0; j < k; ++j) { ret.A[i][j] = 0; for(int h = 0; h < k; ++h) { ret.A[i][j] += A[i][h] * o.A[h][j]; } ret.A[i][j] %= p; } } return ret; } inline int Det2(ll a, ll b, ll c, ll d) { int ret = (a * d - b * c) % p; if (ret < 0) ret += p; return ret; } inline int Det() { if (k == 2) { return Det2(A[0][0], A[0][1], A[1][0], A[1][1]); } else { ll ret = 0; (ret += A[0][0] * A[1][1] * A[2][2]) %= p; (ret += A[0][1] * A[1][2] * A[2][0]) %= p; (ret += A[0][2] * A[1][0] * A[2][1]) %= p; (ret -= A[0][2] * A[1][1] * A[2][0]) %= p; (ret -= A[0][0] * A[1][2] * A[2][1]) %= p; (ret -= A[0][1] * A[1][0] * A[2][2]) %= p; (ret += p) %= p; return ret; } } inline Matrix operator~() { ll d = inv[Det()]; Matrix ret; if (k == 2) { ret.A[0][0] = A[1][1] * d % p; ret.A[1][0] = (p - A[1][0]) * d % p; ret.A[1][1] = A[0][0] * d % p; ret.A[0][1] = (p - A[0][1]) * d % p; } else { ret.A[0][0] = Det2(A[1][1], A[1][2], A[2][1], A[2][2]) * d % p; ret.A[1][0] = Det2(A[1][0], A[1][2], A[2][0], A[2][2]) * d % p; ret.A[2][0] = Det2(A[1][0], A[1][1], A[2][0], A[2][1]) * d % p; ret.A[0][1] = Det2(A[0][1], A[0][2], A[2][1], A[2][2]) * d % p; ret.A[1][1] = Det2(A[0][0], A[0][2], A[2][0], A[2][2]) * d % p; ret.A[2][1] = Det2(A[0][0], A[0][1], A[2][0], A[2][1]) * d % p; ret.A[0][2] = Det2(A[0][1], A[0][2], A[1][1], A[1][2]) * d % p; ret.A[1][2] = Det2(A[0][0], A[0][2], A[1][0], A[1][2]) * d % p; ret.A[2][2] = Det2(A[0][0], A[0][1], A[1][0], A[1][1]) * d % p; for (int i = 0; i < 3; ++i) { for (int j = 0; j < 3; ++j) if ((i + j) & 1) { ret.A[i][j] = p - ret.A[i][j]; } } } return ret; } }; Matrix a[p]; void get_inv(const int x) { inv[1] = 1; for (int i = 2; i < x; ++i) inv[i] = 1ll * (p - p / i) * inv[p % i] % p; } int main() { std::ios::sync_with_stdio(false); std::cin.tie(nullptr); get_inv(p); std::cin >> n >> k >> q; for (int i = 0; i < k; ++i) a[0].A[i][i] = 1; for (int i = 1; i <= n; ++i) { for (int j = 0; j < k; ++j) { for (int h = 0; h < k; ++h) { std::cin >> a[i].A[j][h]; } } a[i] = a[i - 1] * a[i]; } int ans = 0; for (int l, r; q; --q) { std::cin >> l >> r; Matrix mul = (~a[l - 1]) * a[r]; for (int i = 0; i < k; ++i) { for (int j = 0; j < k; ++j) { ans ^= mul.A[i][j]; } } } std::cout << ans << std::endl; } ```