矩阵运算

众所周知，数是可以进行加减乘除的，那矩阵为啥不可以呢？

假设现在我们有两个矩阵 $A$ 和 $B$，矩阵大小分别为 $n \times m$ 和 $x \times y$，矩阵元素对 $mod$ 取模。

const int N = 110, mod = 1e9 + 7; // 矩阵大小和模数

struct matrix {
  int a[N][N], n, m; // 存储矩阵和矩阵大小
} ;

基本运算

矩阵加法

令 $A + B = C$。

要求：$n = x$ 并且 $m = y$。

其实很简单，就是一一对应着加就行，即对于 $1\leqslant i \leqslant n, 1\leqslant j \leqslant m$，$C_{i, j} = A_{i, j} + B_{i, j}$。

所以 $C$ 的大小为 $n \times m$，矩阵减法同。

性质

交换律：明显满足，即 $A + B = B + A$。
结合律：明显满足，即 $A + B + C = A + (B + C)$。

单位矩阵

明显，就是一个大小为 $n \times m$ 的全 $0$ 矩阵。

矩阵加法 Code

时间复杂度：$O(n \times m)$。

  matrix operator + (matrix y) { // 重载运算符 +
    matrix z;
    z.n = n, z.m = m;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= m; j++)
        z.a[i][j] = a[i][j] + y.a[i][j];
    return z;
  }

矩阵乘法

令 $A \times B = C$。

要求：$m = x$。

公式：对于 $1 \leqslant i \leqslant n, 1 \leqslant j \leqslant y$:

\[C_{i, j} = \sum_{1 \leqslant k \leqslant m} A_{i, k} \times B_{k, j} \]

似乎看起来没啥用？因为它是在某些特定方面可以进行时间优化(例如用 $O(8 \log n)$ 的时间复杂度求出斐波那契数列的第 $n$ 项)，所以一般就是把它和其他知识点捆绑起来(例如矩阵乘法优化 dp 等)，单独拎出来出的题目很少(大多都是板子)。

性质

交换律：不满足！
- 若 $n = y$，交换以后还能乘法，但矩阵大小会变为 $m \times x$，与原矩乘后答案不同。
- 当然也有可能 $n \neq y$，这样交换以后连矩乘的基本要求都无法满足。
结合律：满足，设还有一个矩阵 $C$，大小为 $y \times z$，$A \times B \times C = D$，则 $D$ 大小为 $n \times z$，而 $A \times (B \times C) = A \times E(\texttt{大小为} x \times z) = D$。
分配律：明显满足，设还有一个矩阵 $C$，大小为 $y \times z$，即 $A \times (B + C) = A \times B + A \times C$。

单位矩阵

Link

提到矩乘，就不得不提到它的单位矩阵。

单位矩阵定义：若这个矩阵乘任意一个矩阵 $C$，得到的结果都等于 $C$，则这个矩阵为单位矩阵。

观察式子，容易推出：当一个大小为 $n \times n$ 的矩阵除了主对角线(从左上到右下)以外的数为 $1$，其他都为 $0$ 时，它就是一个可以作为所有大小为 $n \times m$ 的矩阵的单位矩阵。

时间复杂度：$O(x \times y)$。

  void Clear (int x, int y, int f) { // 把矩阵大小设为 x * y，矩阵元素都为 f
    n = x, m = y;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= m; j++)
        a[i][j] = f;
  }
  void Unit (int x, int y) { // 构造单位矩阵
    Clear(x, y, 0); // 先清空
    for (int i = 1; i <= n; i++)
      a[i][i] = 1; // 只有主对角线为 1
  }

矩阵乘法 Code

时间复杂度：$O(n \times m \times y)$。

  matrix operator * (matrix y) { // 重载运算符 *
    matrix z;
    z.n = n, z.m = y.m;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= y.m; j++) {
        z.a[i][j] = 0;
        for (int k = 1; k <= m; k++)
          z.a[i][j] = (z.a[i][j] + 1ll * a[i][k] * y.a[k][j] % mod) % mod; // 套公式
      }
    return z;
  }

矩阵除法

令 $\frac{A}{B} = C$。

要求：同矩阵乘法。

由于浮点数不可以取模，数里面除法取模等于乘上除数的逆，所以矩阵除法同矩阵除法，不过是把 $B_{i, j}$ 变成 $B_{i, j}$ 在模 $mod$ 意义下的逆。

矩乘与它的特殊运算和应用

P3390 矩乘快速幂

给定一个矩阵 $A$，大小为 $n \times m$，一个指数 $k$ 和一个模数 $mod$，求 $A^k\ \%\ mod$。

要求：矩阵的大小必须为 $n \times n$。

既然矩乘有结合律，那么就可以考虑使用类似于整数快速幂的方法以 $O(\log k)$ 的方式求出 $A^k$，由于一次矩阵乘法需要 $O(n^3)$，所以总时间复杂度为 $O(n ^ 3 \log k)$。

矩乘快速幂 Code

时间复杂度：$O(n ^ 3 \log k)$。

matrix MATqmod (matrix x, ll y) { // 矩阵乘法快速幂
  matrix sum;
  sum.Unit(x.n, x.m); // 初始化为单位矩阵
  while (y) sum = (y & 1 ? sum * x : sum), x = x * x, y >>= 1;
  return sum;
}

P1349 广义斐波那契数列

P1962 斐波那契数列

广义的斐波那契数列是指形如 $a_n=p\times a_{n-1}+q\times a_{n-2}$ 的数列。

今给定数列的两系数 $p$ 和 $q$，以及数列的最前两项 $a_1$ 和 $ a_2$，另给出两个整数 $n$ 和 $m$，试求数列的第 $n$ 项 $a_n \bmod m$。

令 $F_i$ 表示当前情况下斐波那契数列第 $i$ 项的值。

观察一下，把相邻两项构成一个矩阵，即要从 $f_i = \begin{bmatrix}F_i,F_{i - 1}\end{bmatrix}$ 通过操作变为 $f_{i + 1} = \begin{bmatrix}F_{i + 1},F_i\end{bmatrix} = \begin{bmatrix}p \times F_i + q \times F_{i - 1},F_i\end{bmatrix}$，考虑用矩乘完成。

要能与 $\begin{bmatrix}F_i,F_{i - 1}\end{bmatrix}$ 做矩乘，就需要一个大小为 $2 \times k$ 的矩阵，又因为矩乘之后的矩阵仍然大小为 $1 \times 2$，所以 $k = 2$，即要求一个矩阵 $A = \begin{bmatrix}a,b\\c,d\end{bmatrix}$，使得 $A \times f_i = f_{i + 1}$。

套入公式，就是要求 $\begin{cases}a\times F_i + c \times F_{i - 1} = p \times F_i + q \times F_{i - 1}\\b\times F_i + d \times F_{i - 1}=F_i\end{cases}$，可以推出 $\begin{cases}a=p\\b=1\\c=q\\d=0\end{cases}$，现在要求 $F_n$，初始矩阵为$\begin{bmatrix}a_2,a_1\end{bmatrix}$，每次乘矩阵 $A = \begin{bmatrix}p,1\\q,0\end{bmatrix}$ 后都会使矩阵两项下标加 $1$，特判 $n = 1$ 的情况(答案为 $a_1$)，那么答案就是 $\begin{bmatrix}a_2,a_1\end{bmatrix} \times A^{n - 2}$ 取第一个元素即可，需要矩阵乘法快速幂。

普通的其实就相当于 $a_1=a_2=p=q=1$。

Code

  matrix c, d;
  int n, mod;
  c.Clear(2, 2, 0), d.Clear(1, 2, 0);
  cin >> c.a[1][1] >> c.a[2][1] >> d.a[1][2] >> d.a[1][1] >> n >> mod;
  if (n == 1) {
    cout << d.a[1][2];
    return 0;
  }
  c.a[1][2] = 1, d = d * qmod(c, n - 2);
  cout << d.a[1][1];

完整代码

点击查看最新代码

#include <bits/stdc++.h>
#define _1 (__int128)1

using namespace std;
using ll = long long;

void FileIO (string s) {
  freopen((s + ".in").c_str(), "r", stdin);
  freopen((s + ".out").c_str(), "w", stdout);
}

const int mod = 1e9 + 7, MATsize = 110; // 矩阵大小和模数

struct matrix {
  int a[MATsize][MATsize], n, m; // 存储矩阵和矩阵大小
  void Clear (int x, int y, int f) { // 把矩阵大小设为 x * y，矩阵元素都为 f
    n = x, m = y;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= m; j++)
        a[i][j] = f;
  }
  void Unit (int x, int y) { // 构造单位矩阵
    Clear(x, y, 0); // 先清空
    for (int i = 1; i <= n; i++)
      a[i][i] = 1; // 只有主对角线为 1
  }
  matrix operator + (matrix y) { // 重载运算符 +
    matrix z;
    z.n = n, z.m = m;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= m; j++)
        z.a[i][j] = a[i][j] + y.a[i][j];
    return z;
  }
  matrix operator * (matrix y) { // 重载运算符 *
    matrix z;
    z.n = n, z.m = y.m;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= y.m; j++) {
        z.a[i][j] = 0;
        for (int k = 1; k <= m; k++)
          z.a[i][j] = (z.a[i][j] + 1ll * a[i][k] * y.a[k][j] % mod) % mod; // 套公式
      }
    return z;
  }
} ;

matrix MATqmod (matrix x, ll y) { // 矩阵乘法快速幂
  matrix sum;
  sum.Unit(x.n, x.m); // 初始化为单位矩阵
  while (y) sum = (y & 1 ? sum * x : sum), x = x * x, y >>= 1;
  return sum;
}

signed main () {
  ios::sync_with_stdio(0), cin.tie(0);
  // FileIO("");
  return 0;
}

点击查看无注释代码

#include <bits/stdc++.h>
#define _1 (__int128)1

using namespace std;
using ll = long long;

void FileIO (string s) {
  freopen((s + ".in").c_str(), "r", stdin);
  freopen((s + ".out").c_str(), "w", stdout);
}

const int mod = 1e9 + 7, MATsize = 110;

struct matrix {
  int a[MATsize][MATsize], n, m;
  void Clear (int x, int y, int f) {
    n = x, m = y;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= m; j++)
        a[i][j] = f;
  }
  void Unit (int x, int y) {
    Clear(x, y, 0);
    for (int i = 1; i <= n; i++)
      a[i][i] = 1;
  }
  matrix operator + (matrix y) {
    matrix z;
    z.n = n, z.m = m;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= m; j++)
        z.a[i][j] = a[i][j] + y.a[i][j];
    return z;
  }
  matrix operator * (matrix y) {
    matrix z;
    z.n = n, z.m = y.m;
    for (int i = 1; i <= n; i++)
      for (int j = 1; j <= y.m; j++) {
        z.a[i][j] = 0;
        for (int k = 1; k <= m; k++)
          z.a[i][j] = (z.a[i][j] + 1ll * a[i][k] * y.a[k][j] % mod) % mod;
      }
    return z;
  }
} ;

matrix MATqmod (matrix x, ll y) {
  matrix sum;
  sum.Unit(x.n, x.m);
  while (y) sum = (y & 1 ? sum * x : sum), x = x * x, y >>= 1;
  return sum;
}

signed main () {
  ios::sync_with_stdio(0), cin.tie(0);
  // FileIO("");
  return 0;
}

posted @ 2023-08-05 11:51 wnsyou 阅读(81) 评论(0) 编辑收藏举报

刷新页面返回顶部

wnsyou の blog

矩阵运算

矩阵运算

基本运算

矩阵加法

性质

单位矩阵

矩阵加法 Code

矩阵乘法

性质

单位矩阵

矩阵乘法 Code

矩阵除法

矩乘与它的特殊运算和应用

P3390 矩乘快速幂

矩乘快速幂 Code

P1349 广义斐波那契数列

Code

完整代码

公告