跳转至

矩陣

本文介紹線性代數中一個非常重要的內容——矩陣(Matrix),主要講解矩陣的性質、運算,以及矩陣乘法的一些應用。

向量與矩陣

在線性代數中,向量分為列向量和行向量。

Warning

在中國台灣地區關於「列」與「行」的翻譯,恰好與中國大陸地區相反。在 OI Wiki 按照中國大陸地區的習慣,採用列(column)與行(row)的翻譯。

線性代數的主要研究對象是列向量,約定使用粗體小寫字母表示列向量。在用到大量向量與矩陣的線性代數中,不引起混淆的情況下,在手寫時,字母上方的向量記號可以省略不寫。

向量也是特殊的矩陣。如果想要表示行向量,需要在粗體小寫字母右上方寫轉置記號。行向量在線性代數中一般表示方程。

引入

矩陣的引入來自於線性方程組。與向量類似,矩陣體現了一種對數據「打包處理」的思想。

例如,將線性方程組:

\[ \begin{equation} \begin{cases} 7x_1+8x_2+9x_3=13 \\ 4x_1+5x_2+6x_3=12 \\ x_1+2x_2+3x_3=11 \end{cases} \end{equation} \]

一般用圓括號或方括號表示矩陣。將上述係數抽出來,寫成矩陣乘法的形式:

\[ \begin{equation} \begin{pmatrix} 7 & 8 & 9 \\ 4 & 5 & 6 \\ 1 & 2 & 3 \end{pmatrix}\begin{pmatrix} x_1 \\ x_2 \\ x_3 \end{pmatrix}=\begin{pmatrix} 13 \\ 12 \\ 11 \end{pmatrix} \end{equation} \]

簡記為:

\[ Ax=b \]

即未知數列向量 x,左乘一個矩陣 A,得到列向量 b。這個式子可以認為是線性代數的基本形式。

線性代數主要研究的運算模型是內積。內積是先相乘再相加,是行向量左乘列向量,得到一個數的過程。

矩陣乘法是內積的拓展。矩陣乘法等價於左邊矩陣抽出一行,與右邊矩陣抽出一列進行內積,得到結果矩陣的對應元素,口訣「左行右列」。

當研究對象是右邊的列向量時,矩陣乘法相當於對列向量進行左乘。在左乘的觀點下,矩陣就是對列向量的變換,將矩陣乘法中右邊矩陣的每一個列向量進行變換,對應地得到結果矩陣中每一個列向量。

矩陣可以對一個列向量進行變換,也可以對一組列向量進行「打包」變換,甚至可以對整個空間——即全體列向量進行變換。當矩陣被視為對整個空間變換的時候,也就脱離了空間,成為了純粹變換的存在。

定義

對於矩陣 \(A\),主對角線是指 \(A_{i,i}\) 的元素。

一般用 \(I\) 來表示單位矩陣,就是主對角線上為 1,其餘位置為 0。

同型矩陣

兩個矩陣,行數與列數對應相同,稱為同型矩陣。

方陣

行數等於列數的矩陣稱為方陣。方陣是一種特殊的矩陣。對於「\(n\) 階矩陣」的習慣表述,實際上講的是 \(n\) 階方陣。階數相同的方陣為同型矩陣。

研究方程組、向量組、矩陣的秩的時候,使用一般的矩陣。研究特徵值和特徵向量、二次型的時候,使用方陣。

主對角線

方陣中行數等於列數的元素構成主對角線。

對稱矩陣

如果方陣的元素關於主對角線對稱,即對於任意的 \(i\)\(j\)\(i\)\(j\) 列的元素與 \(j\)\(i\) 列的元素相等,則將方陣稱為對稱矩陣。

對角矩陣

主對角線之外的元素均為 \(0\) 的方陣稱為對角矩陣,一般記作:

\[ \operatorname{diag}\{\lambda_1,\cdots,\lambda_n\} \]

式中的 \(\lambda_1,\cdots,\lambda_n\) 是主對角線上的元素。

對角矩陣是對稱矩陣。

如果對角矩陣的元素均為 \(1\),稱為單位矩陣,記為 \(I\)。只要乘法可以進行,無論形狀,任何矩陣乘單位矩陣仍然保持不變。

三角矩陣

如果方陣主對角線左下方的元素均為 \(0\),稱為上三角矩陣。如果方陣主對角線右上方的元素均為 \(0\),稱為下三角矩陣。

兩個上(下)三角矩陣的乘積仍然是上(下)三角矩陣。如果對角線元素均非 \(0\),則上(下)三角矩陣可逆,逆也是上(下)三角矩陣。

單位三角矩陣

如果上三角矩陣 \(A\) 的對角線全為 \(1\),則稱 \(A\) 是單位上三角矩陣。如果下三角矩陣 \(A\) 的對角線全為 \(1\),則稱 \(A\) 是單位下三角矩陣。

兩個單位上(下)三角矩陣的乘積仍然是單位上(下)三角矩陣,單位上(下)三角矩陣的逆也是單位上(下)三角矩陣。

運算

矩陣的線性運算

矩陣的線性運算分為加減法與數乘,它們均為逐個元素進行。只有同型矩陣之間可以對應相加減。

矩陣的轉置

矩陣的轉置,就是在矩陣的右上角寫上轉置「T」記號,表示將矩陣的行與列互換。

對稱矩陣轉置前後保持不變。

矩陣乘法

矩陣的乘法是向量內積的推廣。

矩陣相乘只有在第一個矩陣的列數和第二個矩陣的行數相同時才有意義。

\(A\)\(P \times M\) 的矩陣,\(B\)\(M \times Q\) 的矩陣,設矩陣 \(C\) 為矩陣 \(A\)\(B\) 的乘積,

其中矩陣 \(C\) 中的第 \(i\) 行第 \(j\) 列元素可以表示為:

\[ C_{i,j} = \sum_{k=1}^MA_{i,k}B_{k,j} \]

在矩陣乘法中,結果 \(C\) 矩陣的第 \(i\) 行第 \(j\) 列的數,就是由矩陣 \(A\)\(i\)\(M\) 個數與矩陣 \(B\)\(j\)\(M\) 個數分別 相乘再相加 得到的。這裏的 相乘再相加,就是向量的內積。乘積矩陣中第 \(i\) 行第 \(j\) 列的數恰好是乘數矩陣 \(A\)\(i\) 個行向量與乘數矩陣 \(B\)\(j\) 個列向量的內積,口訣為 左行右列

線性代數研究的向量多為列向量,根據這樣的對矩陣乘法的定義方法,經常研究對列向量左乘一個矩陣的左乘運算,同時也可以在這裏看出「打包處理」的思想,同時處理很多個向量內積。

矩陣乘法滿足結合律,不滿足一般的交換律。

利用結合律,矩陣乘法可以利用 快速冪 的思想來優化。

在比賽中,由於線性遞推式可以表示成矩陣乘法的形式,也通常用矩陣快速冪來求線性遞推數列的某一項。

優化

首先對於比較小的矩陣,可以考慮直接手動展開循環以減小常數。

可以重新排列循環以提高空間局部性,這樣的優化不會改變矩陣乘法的時間複雜度,但是會在得到常數級別的提升。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
// 以下文的參考代碼為例
mat operator*(const mat& T) const {
  mat res;
  for (int i = 0; i < sz; ++i)
    for (int j = 0; j < sz; ++j)
      for (int k = 0; k < sz; ++k) {
        res.a[i][j] += mul(a[i][k], T.a[k][j]);
        res.a[i][j] %= MOD;
      }
  return res;
}

// 不如
mat operator*(const mat& T) const {
  mat res;
  int r;
  for (int i = 0; i < sz; ++i)
    for (int k = 0; k < sz; ++k) {
      r = a[i][k];
      for (int j = 0; j < sz; ++j)
        res.a[i][j] += T.a[k][j] * r, res.a[i][j] %= MOD;
    }
  return res;
}

方陣的逆

方陣 \(A\) 的逆矩陣 \(P\) 是使得 \(A \times P = I\) 的矩陣。

逆矩陣不一定存在。如果存在,可以使用 高斯消元 進行求解。

方陣的行列式

行列式是方陣的一種運算。

參考代碼

一般來説,可以用一個二維數組來模擬矩陣。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
struct mat {
  LL a[sz][sz];

  mat() { memset(a, 0, sizeof a); }

  mat operator-(const mat& T) const {
    mat res;
    for (int i = 0; i < sz; ++i)
      for (int j = 0; j < sz; ++j) {
        res.a[i][j] = (a[i][j] - T.a[i][j]) % MOD;
      }
    return res;
  }

  mat operator+(const mat& T) const {
    mat res;
    for (int i = 0; i < sz; ++i)
      for (int j = 0; j < sz; ++j) {
        res.a[i][j] = (a[i][j] + T.a[i][j]) % MOD;
      }
    return res;
  }

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

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

看待線性方程組的兩種視角

看待矩陣 A,或者變換 A,有兩種視角。

第一種觀點:按行看,觀察 A 的每一行。這樣一來把 A 看作方程組。於是就有了消元法解方程的過程。

第二種觀點:按列看,觀察 A 的每一列。A 本身也是由列向量構成的。此時相當於把變換 A 本身看成了列向量組,而 x 是未知數係數,思考 A 當中的這組列向量能不能配上未知數,湊出列向量 b。

例如,文章開頭的例子變為:

\[ \begin{equation} \begin{pmatrix} 7 \\ 4 \\ 1 \end{pmatrix}x_1+\begin{pmatrix} 8 \\ 5 \\ 2 \end{pmatrix}x_2+\begin{pmatrix} 9 \\ 6 \\ 3 \end{pmatrix}x_3=\begin{pmatrix} 13 \\ 12 \\ 11 \end{pmatrix} \end{equation} \]

解方程變為研究,是否可以通過調整三個係數 x,使得給定的三個基向量能夠湊出結果的向量。

按列看比按行看更新穎。在按列看的視角下,可以研究線性無關與線性相關。

矩陣乘法的應用

矩陣加速遞推

斐波那契數列(Fibonacci Sequence) 為例。在斐波那契數列當中,\(F_1 = F_2 = 1\)\(F_i = F_{i - 1} + F_{i - 2}(i \geq 3)\)

如果有一道題目讓你求斐波那契數列第 \(n\) 項的值,最簡單的方法莫過於直接遞推了。但是如果 \(n\) 的範圍達到了 \(10^{18}\) 級別,遞推就不行了,此時我們可以考慮矩陣加速遞推。

根據斐波那契數列 遞推公式的矩陣形式:

\[ \begin{bmatrix} F_{n-1} & F_{n-2} \end{bmatrix} \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} = \begin{bmatrix} F_n & F_{n-1} \end{bmatrix} \]

定義初始矩陣 \(\text{ans} = \begin{bmatrix}F_2 & F_1\end{bmatrix} = \begin{bmatrix}1 & 1\end{bmatrix}, \text{base} = \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix}\)。那麼,\(F_n\) 就等於 \(\text{ans} \text{base}^{n-2}\) 這個矩陣的第一行第一列元素,也就是 \(\begin{bmatrix}1 & 1\end{bmatrix} \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix}^{n-2}\) 的第一行第一列元素。

注意

矩陣乘法不滿足交換律,所以一定不能寫成 \(\begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix}^{n-2} \begin{bmatrix}1 & 1\end{bmatrix}\) 的第一行第一列元素。另外,對於 \(n \leq 2\) 的情況,直接輸出 \(1\) 即可,不需要執行矩陣快速冪。

為什麼要乘上 \(\text{base}\) 矩陣的 \(n-2\) 次方而不是 \(n\) 次方呢?因為 \(F_1, F_2\) 是不需要進行矩陣乘法就能求的。也就是説,如果只進行一次乘法,就已經求出 \(F_3\) 了。如果還不是很理解為什麼冪是 \(n-2\),建議手算一下。

下面是求斐波那契數列第 \(n\) 項對 \(10^9+7\) 取模的示例代碼(核心部分)。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
const int mod = 1000000007;

struct Matrix {
  int a[3][3];

  Matrix() { memset(a, 0, sizeof a); }

  Matrix operator*(const Matrix &b) const {
    Matrix res;
    for (int i = 1; i <= 2; ++i)
      for (int j = 1; j <= 2; ++j)
        for (int k = 1; k <= 2; ++k)
          res.a[i][j] = (res.a[i][j] + a[i][k] * b.a[k][j]) % mod;
    return res;
  }
} ans, base;

void init() {
  base.a[1][1] = base.a[1][2] = base.a[2][1] = 1;
  ans.a[1][1] = ans.a[1][2] = 1;
}

void qpow(int b) {
  while (b) {
    if (b & 1) ans = ans * base;
    base = base * base;
    b >>= 1;
  }
}

int main() {
  int n = read();
  if (n <= 2) return puts("1"), 0;
  init();
  qpow(n - 2);
  println(ans.a[1][1] % mod);
}

這是一個稍微複雜一些的例子。

\[ \begin{gathered} f_{1} = f_{2} = 0\\ f_{n} = 7f_{n-1}+6f_{n-2}+5n+4\times 3^n \end{gathered} \]

我們發現,\(f_n\)\(f_{n-1}, f_{n-2}, n\) 有關,於是考慮構造一個矩陣描述狀態。

但是發現如果矩陣僅有這三個元素 \(\begin{bmatrix}f_n& f_{n-1}& n\end{bmatrix}\) 是難以構造出轉移方程的,因為乘方運算和 \(+1\) 無法用矩陣描述。

於是考慮構造一個更大的矩陣。

\[ \begin{bmatrix}f_n& f_{n-1}& n& 3^n & 1\end{bmatrix} \]

我們希望構造一個遞推矩陣可以轉移到

\[ \begin{bmatrix} f_{n+1}& f_{n}& n+1& 3^{n+1} & 1 \end{bmatrix} \]

轉移矩陣即為

\[ \begin{bmatrix} 7 & 1 & 0 & 0 & 0\\ 6 & 0 & 0 & 0 & 0\\ 5 & 0 & 1 & 0 & 0\\ 12 & 0 & 0 & 3 & 0\\ 5 & 0 & 1 & 0 & 1 \end{bmatrix} \]

矩陣表達修改

「THUSCH 2017」大魔法師

大魔法師小 L 製作了 \(n\) 個魔力水晶球,每個水晶球有水、火、土三個屬性的能量值。小 L 把這 \(n\) 個水晶球在地上從前向後排成一行,然後開始今天的魔法表演。

我們用 \(A_i,\ B_i,\ C_i\) 分別表示從前向後第 \(i\) 個水晶球(下標從 \(1\) 開始)的水、火、土的能量值。

小 L 計劃施展 \(m\) 次魔法。每次,他會選擇一個區間 \([l, r]\),然後施展以下 \(3\) 大類、\(7\) 種魔法之一:

  1. 魔力激發:令區間裏每個水晶球中 特定屬性 的能量爆發,從而使另一個 特定屬性 的能量增強。具體來説,有以下三種可能的表現形式:

    • 火元素激發水元素能量:令 \(A_i = A_i + B_i\)
    • 土元素激發火元素能量:令 \(B_i = B_i + C_i\)
    • 水元素激發土元素能量:令 \(C_i = C_i + A_i\)

      需要注意的是,增強一種屬性的能量並不會改變另一種屬性的能量,例如 \(A_i = A_i + B_i\) 並不會使 \(B_i\) 增加或減少。

  2. 魔力增強:小 L 揮舞法杖,消耗自身 \(v\) 點法力值,來改變區間裏每個水晶球的 特定屬性 的能量。具體來説,有以下三種可能的表現形式:

    • 火元素能量定值增強:令 \(A_i = A_i + v\)
    • 水元素能量翻倍增強:令 \(B_i=B_i \cdot v\)
    • 土元素能量吸收融合:令 \(C_i = v\)
  3. 魔力釋放:小 L 將區間裏所有水晶球的能量聚集在一起,融合成一個新的水晶球,然後送給場外觀眾。生成的水晶球每種屬性的能量值等於區間內所有水晶球對應能量值的代數和。需要注意的是,魔力釋放的過程不會真正改變區間內水晶球的能量

值得一提的是,小 L 製造和融合的水晶球的原材料都是定製版的 OI 工廠水晶,所以這些水晶球有一個能量閾值 \(998244353\)。當水晶球中某種屬性的能量值大於等於這個閾值時,能量值會自動對閾值取模,從而避免水晶球爆炸。

小 W 為小 L(唯一的)觀眾,圍觀了整個表演,並且收到了小 L 在表演中融合的每個水晶球。小 W 想知道,這些水晶球藴涵的三種屬性的能量值分別是多少。

由於矩陣的結合律和分配律成立,單點修改可以自然地推廣到區間,即推出矩陣後直接用線段樹維護區間矩陣乘積即可。

下面將舉幾個例子。

\(A_i = A_i + v\) 的轉移

\[ \begin{bmatrix} A & B & C & 1 \end{bmatrix} \begin{bmatrix} 1 & 0 & 0 & 0\\ 0 & 1 & 0 & 0\\ 0 & 0 & 1 & 0\\ v & 0 & 0 & 1\\ \end{bmatrix}= \begin{bmatrix} A+v & B & C & 1\\ \end{bmatrix} \]

\(B_i=B_i \cdot v\) 的轉移

\[ \begin{bmatrix} A & B & C & 1 \end{bmatrix} \begin{bmatrix} 1 & 0 & 0 & 0\\ 0 & v & 0 & 0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1\\ \end{bmatrix}= \begin{bmatrix} A & B \cdot v & C & 1\\ \end{bmatrix} \]

「LibreOJ 6208」樹上詢問

有一棵 \(n\) 節點的樹,根為 \(1\) 號節點。每個節點有兩個權值 \(k_i, t_i\),初始值均為 \(0\)

給出三種操作:

  1. \(\operatorname{Add}( x , d )\) 操作:將 \(x\) 到根的路徑上所有點的 \(k_i\leftarrow k_i + d\)
  2. \(\operatorname{Mul}( x , d )\) 操作:將 \(x\) 到根的路徑上所有點的 \(t_i\leftarrow t_i + d \times k_i\)
  3. \(\operatorname{Query}( x )\) 操作:詢問點 \(x\) 的權值 \(t_x\)

    \(n,~m \leq 100000, ~-10 \leq d \leq 10\)

若直接思考,下放操作和維護信息並不是很好想。但是矩陣可以輕鬆地表達。

\[ \begin{aligned} \begin{bmatrix}k & t & 1 \end{bmatrix} \begin{bmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ d & 0 & 1 \end{bmatrix} &= \begin{bmatrix}k+d & t & 1 \end{bmatrix}\\ \begin{bmatrix}k & t & 1 \end{bmatrix} \begin{bmatrix} 1 & d & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{bmatrix} &= \begin{bmatrix}k & t+d \times k & 1 \end{bmatrix} \end{aligned} \]

定長路徑統計

問題描述

給一個 \(n\) 階有向圖,每條邊的邊權均為 \(1\),然後給一個整數 \(k\),你的任務是對於所有點對 \((u,v)\) 求出從 \(u\)\(v\) 長度為 \(k\) 的路徑的數量(不一定是簡單路徑,即路徑上的點或者邊可能走多次)。

我們將這個圖用鄰接矩陣 \(G\)(對於圖中的邊 \((u\to v)\),令 \(G[u,v]=1\),其餘為 \(0\) 的矩陣;如果有重邊,則設 \(G[u,v]\) 為重邊的數量)表示這個有向圖。下述算法同樣適用於圖有自環的情況。

顯然,該鄰接矩陣對應 \(k=1\) 時的答案。

假設我們知道長度為 \(k\) 的路徑條數構成的矩陣,記為矩陣 \(C_k\),我們想求 \(C_{k+1}\)。顯然有 DP 轉移方程

\[ C_{k+1}[i,j] = \sum_{p = 1}^{n} C_k[i,p] \cdot G[p,j] \]

我們可以把它看作矩陣乘法的運算,於是上述轉移可以描述為

\[ C_{k+1} = C_k \cdot G \]

那麼把這個遞推式展開可以得到

\[ C_k = \underbrace{G \cdot G \cdots G}_{k \text{ 次}} = G^k \]

要計算這個矩陣冪,我們可以使用快速冪(二進制取冪)的思想,在 \(O(n^3 \log k)\) 的複雜度內計算結果。

定長最短路

問題描述

給你一個 \(n\) 階加權有向圖和一個整數 \(k\)。對於每個點對 \((u,v)\) 找到從 \(u\)\(v\) 的恰好包含 \(k\) 條邊的最短路的長度。(不一定是簡單路徑,即路徑上的點或者邊可能走多次)

我們仍構造這個圖的鄰接矩陣 \(G\)\(G[i,j]\) 表示從 \(i\)\(j\) 的邊權。如果 \(i,j\) 兩點之間沒有邊,那麼 \(G[i,j]=\infty\)。(有重邊的情況取邊權的最小值)

顯然上述矩陣對應 \(k=1\) 時問題的答案。我們仍假設我們知道 \(k\) 的答案,記為矩陣 \(L_k\)。現在我們想求 \(k+1\) 的答案。顯然有轉移方程

\[ L_{k+1}[i,j] = \min_{1\le p \le n} \left\{L_k[i,p] + G[p,j]\right\} \]

事實上我們可以類比矩陣乘法,你發現上述轉移只是把矩陣乘法的乘積求和變成相加取最小值,於是我們定義這個運算為 \(\odot\),即

\[ A \odot B = C~~\Longleftrightarrow~~C[i,j]=\min_{1\le p \le n}\left\{A[i,p] + B[p,j]\right\} \]

於是得到

\[ L_{k+1} = L_k \odot G \]

展開遞推式得到

\[ L_k = \underbrace{G \odot \ldots \odot G}_{k\text{ 次}} = G^{\odot k} \]

我們仍然可以用矩陣快速冪的方法計算上式,因為它顯然是具有結合律的。時間複雜度 \(O(n^3 \log k)\)

限長路徑計數/最短路

上述算法只適用於邊數固定的情況。然而我們可以改進算法以解決邊數小於等於 \(k\) 的情況。具體地,考慮以下問題:

問題描述

給一個 \(n\) 階有向圖,邊權為 \(1\),然後給一個整數 \(k\),你的任務是對於每個點對 \((u,v)\) 找到從 \(u\)\(v\) 長度小於等於 \(k\) 的路徑的數量(不一定是簡單路徑,即路徑上的點或者邊可能走多次)。

我們簡單修改一下這個圖,我們給每一個結點加一個權值為 \(1\) 的自環。這樣走的時侯就可以走自環,相當於原地走。這樣就包含了小於等於 \(k\) 的情況。修改後再做矩陣快速冪即可。(即使這個圖在修改之前就有自環,該算法仍是成立的)。

同樣的方法可以用於求邊數小於等於 \(k\) 的最短路,即加一個邊權為 \(0\) 的自環。

習題