自动搬运

查看原文

来自洛谷，原作者为

	一扶苏一 休息结束。邮箱 yifusuyi@qq.com

搬运于2025-08-24 21:14:08，当前版本为作者最后更新于2018-06-26 18:45:57，作者可能在搬运后再次修改，您可在原文处查看最新版

自动搬运只会搬运当前题目点赞数最高的题解，您可前往洛谷题解查看更多

以下是正文

数列前缀和 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$ 只有 $2$ 和 $3$，所以也可以手搓一下伴随矩阵，下面介绍一下伴随矩阵法求逆矩阵。

对 $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$ 的伴随矩阵 $A^*$ 也是一个 $n$ 阶矩阵，其第 $i$ 行第 $j$ 列为 $a_{j, i}$ 的代数余子式 $A_{j, i}$（注意这里下标是反着的，即伴随矩阵的 $i$ 行 $j$ 列是原矩阵 $j$ 行 $i$ 列的代数余子式）。可以证明，$A^{-1} = \frac{A^*}{|A|}$，其中 $|A|$ 是 $A$ 的行列式。

对于二阶矩阵 $A = \begin{pmatrix} a & b \\ c & d\end{pmatrix}$，直接套用上述结论得到 $A^{-1} = \frac{1}{ad - bc}\begin{pmatrix} d & -b \\ -c & a\end{pmatrix}$；对于三阶矩阵，其行列式只有六项，所有的代数余子式均为二阶行列式，都可以轻松算出。

#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;
}