acm模板-线性代数

发布于 2021-10-20  1661 次阅读


线性代数

矩阵操作 & 高斯-约旦消元法 & 矩阵求逆

struct mat {
  vector<vector<LL>> a;

  inline mat(int n, int m) {
    a.resize(n);
    for (int i = 0; i < n; ++i) a[i].resize(m, 0);
  }

  vector<LL> &operator[](int x) { return a[x]; }

  inline mat operator*(mat T) const {
    mat res(a.size(), T[0].size());
    for (int i = 0; i < a.size(); ++i)
      for (int k = 0; k < a[0].size(); ++k) {
        auto r = a[i][k];
        for (int j = 0; j < T[0].size(); ++j)
          res[i][j] = (res[i][j] + T[k][j] * r) % MOD;
      }
    return res;
  }

  inline mat operator^(LL x) const {
    mat res(a.size(), a[0].size()), bas = *this;
    for (int i = 0; i < a.size(); ++i) res[i][i] = 1;
    while (x) {
      if (x & 1) res = res * bas;
      bas = bas * bas;
      x >>= 1;
    }
    return res;
  }

  bool gauss() {
    for (int i = 0; i < a.size(); ++i) {
      auto mx = i;
      for (int j = i + 1; j < a.size(); ++j)
        if (abs(a[j][i]) > abs(a[mx][i])) mx = j;
      swap(a[i], a[mx]);
      if (abs(a[i][i]) < eps) return false;
      for (int j = 0; j < a.size(); ++j)
        if (j != i) {
          auto tmp = a[j][i] / a[i][i];
          for (int k = i; k < a[i].size(); ++k)
            a[j][k] = a[j][k] - a[i][k] * tmp;
        }
    }
    return true;
  }

  bool inverse() {
    mat res(a.size(), a[0].size() * 2);
    for (int i = 0; i < a.size(); ++i) {
      for (int j = 0; j < a[0].size(); ++j) res[i][j] = a[i][j];
      res[i][i + a[0].size()] = 1;
    }
    if (!res.gauss()) return false;
    for (int i = 0; i < a.size(); ++i) {
      auto tmp = ksm(res[i][i], MOD - 2, MOD);
      for (int j = 0; j < a[i].size(); ++j)
        a[i][j] = (res[i][j + a[i].size()] * tmp) % MOD;
    }
    return true;
  }
};

高斯消元法解异或方程组

bitset<1010> matrix[2010];  // matrix[1~n]:增广矩阵,0 位置为常数
vector<bool> GaussElimination(
    int n, int m)  // n 为未知数个数,m 为方程个数,返回方程组的解(多解 /
                   // 无解返回一个空的 vector)
{
  for (int i = 1; i <= n; i++) {
    int cur = i;
    while (cur <= m && !matrix[cur].test(i)) cur++;
    if (cur > m) return vector<bool>(0);
    if (cur != i) swap(matrix[cur], matrix[i]);
    for (int j = 1; j <= m; j++)
      if (i != j && matrix[j].test(i)) matrix[j] ^= matrix[i];
  }
  vector<bool> ans(n + 1, 0);
  for (int i = 1; i <= n; i++) ans[i] = matrix[i].test(0);
  return ans;
}
Hello, world!
最后更新于 2021-10-20