這是校內講課課件。其實這些程式碼我覺得並沒有必要放上，不過作為許多題解的合集，海蝕放上了。

矩陣簡介

相信大家都會，酌情跳過。

定義

矩陣就是一個填著若干個數的表格。不是什麼高深東西。。。

大概長這樣：

\[\begin{pmatrix} 1&2&3&5\\ 2&1&3&1\\ 1&3&3&2\\ \end{pmatrix} \]

運算

加法運算

兩個大小相等的矩陣才能相加。對位相加即可。計算複雜度為 \(\mathcal O(n\times m)\) 。

乘法運算

你知道矩陣嗎？你知道乘法嗎？那你不就知道矩陣乘法了嘛。

兩個矩陣可以相乘，當且僅當第一個矩陣的列數與第二個矩陣的行數相等。

大小分別為 \(n\times m\) 和 \(m\times k\) 的矩陣 \(a\) 和 \(b\) 的乘法這樣定義（相乘的結果是 \(c\) ，是一個 \(n\times k\) 的矩陣）：

\[c_{i,j}=\sum_{k=1}^{m}a_{i,k}\cdot b_{k,j} \]

也就是說，第一個矩陣的第一行會和第二個矩陣的第一列對位相乘，然後再加起來得到 \(c_{1,1}\) 。

舉個例子：

\[\begin{pmatrix} a_{1,1}&a_{1,2} \end{pmatrix} \times \begin{pmatrix} b_{1,1}\\ b_{2,1}\\ \end{pmatrix} = \begin{pmatrix} a_{1,1}\cdot b_{1,1}+a_{1,2}\cdot b_{2,1} \end{pmatrix} \]

計算時間複雜度為 \(\mathcal O(n\times m\times k)\) 。

運算規律

加法

加法顯然符合交換律、結合律。

乘法

乘法不符合交換律，符合結合律。

證明

上面兩式子都等於

“矩陣運算結合律”本質上就是“獨立指標求和順序無關”。

——來自知乎使用者 瓦拉莫古力斯 的解答

矩陣快速冪簡介

引入

求斐波那契數列的第 \(n\) 項，\(n\le 10^{18}\) 。模 998244353。

解法

你發現這和矩陣快速冪沒有任何關係，因此你暴斃了。

設 \(f_i\) 表示斐波那契數的 \(i\) 位是多少。

首先可以發現

\[\begin{pmatrix} f_i& f_{i+1} \end{pmatrix} \cdot \begin{pmatrix} 0&1\\ 1&1\\ \end{pmatrix} = \begin{pmatrix} f_{i+1}&f_{i+2} \end{pmatrix} \]

所以

\[\begin{pmatrix} f_1&f_2 \end{pmatrix} \cdot \begin{pmatrix} 0&1\\ 1&1 \end{pmatrix}^{n} = \begin{pmatrix} f_{n+1}&f_{n+2} \end{pmatrix} \]

矩陣快速冪即可。

寫法

這裡介紹一下比較好的矩乘寫法。

struct Mat {
    int A[N][N];
    Mat operator*(const Mat &d) const {
        //寫一個矩陣乘法
    }
    int* operator[](const int i) const { return A[i]; }
    //這樣可以讓mt的第i行j列這樣訪問：mt[i][j]，而不用mt.A[i][j]
    //當然這個可能比較玄妙，可能會出現意外錯誤，雖然我沒有遇到過（但我也只用過一次這麼寫的。。）
}

例題

由於這個知識點比較簡單所以需要放一大堆例題拖延時間。

聖人的數列

題目連結

點選開啟連結

題解

我記得 l 學長去年講過這題

其實這個就是指數上的遞推，就像是 \(f_i=f_{i-1}\cdot f_{i-2}^2\) 這種。只要記錄 \(f_i\) 是 \(f_1\) 的多少次放就好了。

記錄每一項 \(f_1\) 的次數，\(f_2\)的次數，\(f_3\)的次數與 \(c\) 的次數，這樣使用四個矩陣轉移出最終那一項的 \(f_1\) ，\(f_2\) ，\(f_3\) 和 \(c\) 的次數，然後快速冪回去就行了。

總結

• 乘法的本質實際上是指數相加

程式碼

然而並沒有程式碼。講題人並沒有寫...

「SNOI2017」禮物

題目連結

點選開啟連結

題解

其實就是求，\(f_i=2\cdot f_{i-1}+i^k\) 的第 \(n\) 項。

下面是奇技淫巧：

\[\begin{pmatrix} f_{i-1}&i^0&i^1&i^2 \end{pmatrix} \cdot \begin{pmatrix} 2&0&0&0\\ 0&1&1&1\\ 0&0&1&2\\ 1&0&0&1\\ \end{pmatrix} = \begin{pmatrix} f_{i}&(i+1)^0&(i+1)^1&(i+2)^2 \end{pmatrix} \]

其中轉移矩陣的右下角是個楊輝三角。

總結

• 轉移矩陣上放楊輝三角的小 trick

程式碼

#include <cstdio>
#include <cstring>
#define int long long

using namespace std;

const int mod = 1e9 + 7;

struct matrix
{
	int row, column, num[15][15];
	
	matrix()
	{
		row = 0;
		column = 0;
		memset(num, 0, sizeof(num));
	}
	
	void unit(int k)
	{
		row = k; column = k;
		for (int i = 1; i <= k; ++i)
			num[i][i] = 1;
	}
	
	matrix operator * (const matrix &d)
	{
		matrix ret;
		ret.row = row;
		ret.column = d.column;
		for (int i = 1; i <= ret.row; ++i)
			for (int j = 1; j <= ret.column; ++j)
				for (int k = 1; k <= column; ++k)
					ret.num[i][j] = (ret.num[i][j] + num[i][k] * d.num[k][j]) % mod;
		return ret;
	}
} ;

matrix power(matrix base, int ciShu)
{
	matrix ret;
	ret.unit(base.row);
	while (ciShu)
	{
		if (ciShu & 1) ret = ret * base;
		base = base * base;
		ciShu >>= 1;
	}
	return ret;
}

int c(int x, int y)
{
	int ret = 1;
	for (int i = x; i >= x - y + 1; --i)
		ret *= i;
	for (int i = 1; i <= y; ++i)
		ret /= i;
	return ret;
}

void solve()
{
	int n, k;
	scanf("%lld%lld", &n, &k);
	if (n == 1)
	{
		printf("1\n");
		return;
	}
	matrix base;
	base.row = base.column = k + 2;
	//注意這裡變數名弄重過一次！本來寫的是base.row = base.column = k + 2;
	base.num[1][1] = 2;
	for (int i = 2; i <= k + 2; ++i)
		base.num[i][1] = c(k, i - 2);
	base.num[2][2] = 1;
	for (int i = 3; i <= k + 2; ++i)
		for (int j = 2; j <= i; ++j)
			base.num[j][i] = base.num[j][i - 1] + base.num[j - 1][i - 1];
	matrix sn, sn_1;
	sn.row = 1;
	sn.column = k + 2;
	for (int i = 1; i <= k + 2; ++i)
		sn.num[1][i] = 1;
	sn_1 = sn;
	sn = sn * power(base, n - 1);
	sn_1 = sn_1 * power(base, n - 2);
	printf("%lld\n", (((sn.num[1][1] - sn_1.num[1][1]) % mod) + mod) % mod);
}

signed main()
{
	solve();
	return 0;
}

【BNDSOJ1396】有向圖的路徑方案數

題目連結

點選開啟連結

題解

很明顯，這可以dp。

設 \(f_{i,j}\) 表示從 A 出發走 \(i\) 步到達 \(j\) 的方案數。轉移即可。然而這樣並不能過。

由於這是一個遞推式，所以可以使用矩陣乘法優化。發現轉移矩陣是鄰接矩陣。

總結

• 這題本身就是個經典套路。。

程式碼

#include <cstdio>
#include <cstring>
#include <iostream>
#define int long long
#define reg register

using namespace std;

inline char nextchar() {
  static char buf[100000], *p1 = buf, *p2 = buf;
  return (p1 == p2) &&
    (p2 = (p1 = buf) + fread(buf, 1, 100000, stdin), p1 == p2) ? EOF : *p1++;
}

inline int read() {
  int x = 0;
  char ch = 0;
  bool sign = false;
  while (!isdigit(ch)) { sign |= (ch == '-'); ch = nextchar(); }
  while (isdigit(ch)) { x = x * 10 + (ch ^ 48); ch = nextchar(); }
  x = sign ? -x : x;
  return x;
}

const int N = 25, MOD = 1000;
struct matrix
{
	int row, column, num[N][N];
} linjie;
int n, m;

void init(matrix &p)
{
	p.row = 0, p.column = 0;
	memset(p.num, 0, sizeof(p.num));
}

void unit(matrix &p, int wyy)
{
	init(p);
	p.column = wyy, p.row = wyy;
	for (reg int i = 1; i <= n; ++i)
		p.num[i][i] = 1;
}

matrix mul(matrix a, matrix d)
{
	matrix ret;
	init(ret);
	ret.row = a.row;
	ret.column = d.column;
	for (reg int i = 1; i <= a.row; ++i)
		for (reg int j = 1; j <= d.column; ++j) 
			for (reg int k = 1; k <= a.column; ++k)
				ret.num[i][j] += a.num[i][k] * d.num[k][j];
	for (int i = 1; i <= a.row; ++i)
		for (int j = 1; j <= d.column; ++j)
			ret.num[i][j] %= MOD;
	return ret;
}

matrix power(matrix a, int b)
{
	matrix ret;
	init(ret);
	unit(ret, a.row);
	while (b)
	{
		if (b & 1) ret = mul(ret, a);
		a = mul(a, a);
		b >>= 1;
	}
	return ret;
}

inline void solve()
{
	n = read(), m = read();
	init(linjie);
	linjie.row = n; linjie.column = n;
	for (reg int i = 1; i <= m; ++i)
	{
		int u, v;
		u = read(), v = read();
		linjie.num[u][v] = 1;
	}
	int q;
	q = read();
	for (reg int i = 1; i <= q; ++i)
	{
		int a, b, k;
		a = read(), b = read(), k = read();
		matrix ans;
		init(ans);
		ans.row = 1;
		ans.column = n;
		ans.num[1][a] = 1;
		ans = mul(ans, power(linjie, k));
		printf("%lld\n", ans.num[1][b] % MOD);
	}
}

int T;

signed main()
{
	T = read();
	while (T--) solve();
	return 0;
}

【luoguP2151】【SDOI2009】HH去散步

題目連結

點選開啟連結

題解

這又是一個小trick嘍。

考慮邊當點，點當邊。在由點組成的圖中的一條路徑，比如是 \(1\rightarrow 2\rightarrow 3\)，那麼在由邊當成點組成的路徑是這樣的：\(12\rightarrow 23\)（這裡以邊的兩個端點的編號連起來作為一個邊的編號，雖然題中說由重複邊，但這無傷大雅）。那麼考慮一條在原圖中走出來又立刻走回去的路徑，比如 \(1\rightarrow 2\rightarrow 1\)，在邊當成的點組成的圖中路徑是這樣的：\(12\rightarrow 21\)，也就是在同一個結點停留了。所以只要在由邊當成點組成的圖中路徑找長度為 \(t-1\) 的路徑即可。這個很明顯可以鄰接矩陣快速冪做。

總結

• 邊當點，點當邊。有的時候可以通過資料範圍觀察出這樣做。
• 可以用邊當點，點當邊解決“不走上回走來的邊”這種問題。

程式碼

然而並沒有。

【BZOJ4386】【POI2015】Wycieczki

題目連結

點選開啟連結

題目解法

這道題求第 \(k\) 短路。所以考慮使用 A* 。考慮統計出長度 \(\le l\) 的路徑數量，然後二分 \(l\) 。

由於有 1，2，3 三種長度，但是我們只會求邊長為 1 的圖的路徑數量。所以考慮將 \(i\) 拆為 3 個點 \(i,0\) 和 \(i,1\) 和 \(i,2\)，考慮將 \(i\) 到 \(j\) 長度為 k 的路徑從 \(i,k-1\) 連線到 \(j,0\) 。然後再把 \(i,0\) 連向 \(i,1\) ，\(i,1\) 連向 \(i,2\) 。這樣的圖就與原圖等價了。

但是即使這樣我們只會求圖上路徑長度恰好為 \(k\) 的路徑數量。怎麼辦呢？考慮新增一個 0 號點，將所有點向 0 連邊，並在 0 號點上加一個自環。這樣的話所有長度不夠的路徑可以走到 0 號點然後不斷走那個自環。這樣的話時間複雜度 \(\mathcal O(n^3\log^2k)\) 。但是這很明顯過不了。

所以考慮優化。容易發現二分+快速冪可以被倍增代替。具體來講，首先可以 \(n^3\log n\) 求出矩陣的若干次冪，然後每次計算其他所有點到 0 號點的路徑數量數量和，這也就是這張圖上的所有 \(\le k\) 的路徑數量。

總結

• 拆點小 trick
• 統計長度 \(\le k\) 的路徑數量的方法
• 二分+快速冪可以被倍增代替

程式碼

並沒有調好呢。

【CodeForces575A】Fibonotci

題目連結

點選開啟連結

題目解法

考慮到這是個遞推式，因此採用矩陣優化地推。那麼構造以下算式：

\[\begin{pmatrix}f_{0}&f_{1}\end{pmatrix}\cdot \begin{pmatrix}1&s_0\\0&s_1\\\end{pmatrix}\cdots\begin{pmatrix}1&s_n\\0&s_{n+1}\end{pmatrix}=\begin{pmatrix}f_{n+1}& f_{n+2}\end{pmatrix} \]

所以我們只需要求一大堆矩陣的乘積。

由於這是迴圈的，所以很多地方可以直接當成一個迴圈節的矩陣積的若干次方。有時還要乘上字首積與字尾積。不過這樣你會發現你還需要乘上區間積：如果兩個修改在同一個區間內，那麼這兩個修改之間的矩陣乘積還需要計算一下，但是你發現使用字首積計算是不成立的，因為不一定每個字首積都可以求逆。

所以你需要一個線段樹維護區間乘積。時間複雜度 \(\mathcal O(m(\log k+\log n))\) 。不過這種方法我沒有寫過，而且也沒看到任何題解這麼寫的，因此我懷疑有一些我沒考慮到的問題，各位可以嘗試這樣寫一下。

不過還有另外一種方法（這種方法可能有一大堆細節問題需要考慮）：

不帶修改的週期可以直接快速冪。考慮用線段樹維護有修改的週期，每次線上段樹上單點修改一個矩陣，然後再求整體乘積。

下面程式碼是這種方法。

總結

• 字首積無法維護矩陣區間積。但由於資訊可以合併，所以可以使用線段樹維護

程式碼

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;
typedef long long ll;
const ll NM = 5e4 + 5;
typedef long long ll;
ll k, mod, n, m;
ll ml(ll x, ll y) { return 1ll * x * y % mod; }
ll ad(ll x, ll y) { return (x + y > mod) ? (x + y - mod) : (x + y); }
ll dc(ll x, ll y) { return (x - y < 0) ? (x - y + mod) : (x - y); }
struct Mat {
	ll A[2][2];
	Mat() {}
	ll * operator [] (ll i) { return A[i]; }
	const ll * operator [] (ll i) const { return A[i]; }
	void Init() {
		A[0][0] = 1;
		A[0][1] = 0;
		A[1][0] = 0;
		A[1][1] = 1;
	}
	void Clear() { memset(A, 0, sizeof(A)); }
	void print() {
		for (ll i = 0; i < 2; ++i) {
			printf("	");
			for (ll j = 0; j < 2; ++j)
				printf("%lld ", A[i][j]);
			printf("\n");
		}
	}
};
struct LYK {
	ll i, x, y, v;
	LYK() {}
	LYK(ll a, ll b, ll c, ll d) : i(a), x(b), y(c), v(d) {}
} szj[NM << 1];
struct SZJ {
	ll i;
	Mat A;
} chg[NM << 1];
bool cmp(const LYK &x, const LYK &y) { return x.i < y.i; }
Mat operator * (const Mat &x, const Mat &y) {
	Mat ret;
	ret.Clear();
	for (ll i = 0; i < 2; ++i)
		for (ll j = 0; j < 2; ++j)
			for (ll k = 0; k < 2; ++k)
				ret[i][j] = ad(ret[i][j], ml(x[i][k], y[k][j]));
	return ret;
}
ll s[NM], ct, mt;
ll blk[NM << 1];
Mat ltg[NM], bas;
void ksm(Mat &x, Mat y, ll c) {
	if (!c) return ;
	for (; c; c >>= 1, y = y * y)
		if (c & 1) x = x * y;
}
Mat tr[NM << 2];
void Upd(ll o) { tr[o] = tr[o << 1] * tr[o << 1 | 1]; }
void Build(ll o, ll sl, ll sr) {
	if (sl == sr) {
		tr[o] = ltg[sl];
		return ;
	}
	ll mid = (sl + sr) >> 1;
	Build(o << 1, sl, mid);
	Build(o << 1 | 1, mid + 1, sr);
	Upd(o);
}
void Modify(ll o, ll sl, ll sr, ll ps, Mat v) {
	if (sl == sr) {
		tr[o] = v;
		return ;
	}
	ll m = (sl + sr) >> 1;
	if (ps <= m) Modify(o << 1, sl, m, ps, v);
	else Modify(o << 1 | 1, m + 1, sr, ps, v);
	Upd(o);
}
int main() {
	scanf("%lld%lld", &k, &mod);
	if (k == 0) printf("0\n");
	if (k == 1) printf("%d\n", mod != 1);
	if (k == 0 || k == 1) return 0;
	--k;
	scanf("%lld", &n);
	for (ll i = 0; i < n; ++i) {
		scanf("%lld", &s[i]);
		s[i] = s[i] % mod;
	}
	for (ll i = 0; i < n; ++i) {
		ltg[i][0][0] = 0;
		ltg[i][0][1] = s[i];
		ltg[i][1][0] = 1;
		ltg[i][1][1] = s[(i == n - 1) ? 0 : (i + 1)];
	}
	bas.Init(); //base的初始化 
	for (ll i = 0; i < n; ++i)
		bas = bas * ltg[i];
	scanf("%lld", &m);
	for (ll i = 1; i <= m; ++i) {
		ll j, v;
		scanf("%lld%lld", &j, &v);
		szj[++ct] = LYK(j, 0, 1, v);
		szj[++ct] = LYK(j - 1, 1, 1, v);
	}
	sort(szj + 1, szj + ct + 1, cmp);
	for (ll i = 1; i <= ct; ++i) {
		if (i != 1 && szj[i].i == szj[i - 1].i)
			chg[mt].A[szj[i].x][szj[i].y] = szj[i].v;
		else {
			chg[++mt].A = ltg[szj[i].i % n];
			chg[mt].A[szj[i].x][szj[i].y] = szj[i].v;
			chg[mt].i = szj[i].i;
		}
	}
	for (ll i = 1; i <= mt; ++i) blk[i] = chg[i].i / n;
	Mat ans;
	ans.Clear();
	ans[0][0] = 0;
	ans[0][1] = 1;
	ll ed = -1;
	Build(1, 0, n - 1);
	ll i, lst = -1, finb = k / n;
	for (i = 1; i <= mt; ) {
		if (blk[i] >= finb)
			break ;
		ksm(ans, bas, blk[i] - lst - 1);
		ll j = i;
		for (; blk[i] == blk[j]; ++i)
			Modify(1, 0, n - 1, chg[i].i % n, chg[i].A);
		ans = ans * tr[1];
		for (ll tt = j; tt < i; ++tt)
			Modify(1, 0, n - 1, chg[tt].i % n, ltg[chg[tt].i % n]);
		lst = blk[i - 1];
		ed = (lst + 1) * n - 1;
	}
	ll nbblk = (k - ed) / n;
	if (k % n == n - 1)
		--nbblk;
	ksm(ans, bas, nbblk);
	ed = ed + nbblk * n;
	ll res = k - ed;
	ll pt = i;
	for (i = 1; i <= res; ++i) {
		if (pt <= mt && chg[pt].i == ed + i)
			ans = ans * chg[pt++].A;
		else ans = ans * ltg[i - 1];
	}
	printf("%lld\n", ans[0][0]);
	return 0;
}

【CodeForces576D】Flights for Regular Customers

題目連結

點選開啟連結

題目解法

考慮列舉最優狀態是哪些邊解鎖了，哪些邊沒有。這樣真的有 \(2^n\) 種情況嗎？並不是的。將所有邊按照解鎖需要走過邊的數量排序。解鎖邊的順序必然是排好序的這樣。所以考慮只解鎖前 \(i\) 條邊，走到 \(n\) 的最短距離是多少。

所以我們可以這樣：每次列舉新解鎖了哪一條邊，如果當前邊需要走 \(i\) 步才能解鎖，那麼看走 \(i\) 條邊並解鎖這個邊後可以到哪些點，用那些點跑一邊多源 dfs 尋找到終點的最短路，更新答案即可。

然而這時間複雜度 \(\mathcal O(n^3m\log n)\) 並不能過去。瓶頸是矩陣快速冪。由於矩陣是0/1 矩陣，所以使用bitset優化即可。比較無腦的bitset可以維護矩陣的行和列，維護兩個bitset陣列，這樣比較方便寫。當然也可以只維護一個bitset，這樣難理解，而且非常非常沒必要。

時間複雜度為 \(\mathcal O(\frac{n^3m\log n}{\omega})\) 。

總結

• bitset 優化0/1矩陣的乘法。
• 解鎖順序可以去掉一些邊的解鎖狀態。

程式碼

#include <cstdio>
#include <bitset>
#include <queue>
#include <algorithm>
using namespace std;
typedef long long ll;
const int NM = 155;
const ll infll = 0x3f3f3f3f3f3f3f3fll;
typedef bitset<150> bt;
int n, m;
struct Edge {
	int u, v, w;
	inline bool operator < (const Edge &d) const { return w < d.w; }
} eg[NM];
struct Mat {
	bt A[150];
	inline bt & operator [] (const int i) { return A[i]; }
	inline const bt & operator [] (const int i) const { return A[i]; }
	//因為 Mat * 不能儲存 const Mat * 的東西
	//所以必須有一種返回 const Mat * 的呼叫方法 
} A;
bt operator * (const Mat &A, const bt &B) {
	static bt ret;
	for (int i = 0; i < n; ++i) ret[i] = (A[i] & B).any();
	return ret;
}
Mat operator * (const Mat &A, const Mat &B) {
	static Mat ret;
	for (int i = 0; i < n; ++i) ret[i].reset();
	for (int i = 0; i < n; ++i)
		for (int j = 0; j < n; ++j)
			if (A[i][j]) ret[i] |= B[j];
	return ret;
}
void ksm(bt &A, Mat B, int c) {
	for (; c; c >>= 1, B = B * B)
		if (c & 1) A = B * A;
}
bt go;
queue<int> Q;
ll dis[NM];
int main() {
	scanf("%d%d", &n, &m);
	for (int i = 1; i <= m; ++i) {
		scanf("%d%d%d", &eg[i].u, &eg[i].v, &eg[i].w);
		--eg[i].u, --eg[i].v;
	}
	sort(eg + 1, eg + m + 1);
	go[0] = 1;
	ll ans = infll;
	for (int i = 1, t = 0; i <= m; ++i) {
		if (eg[i].w >= ans) break ;
		if (eg[i].w != t) ksm(go, A, eg[i].w - t);
		t = eg[i].w;
		A[eg[i].v][eg[i].u] = 1;
		for (int j = 0; j < n; ++j)
			if (go[j]) Q.push(j), dis[j] = 0;
			else dis[j] = infll;
		while (Q.size()) {
			int u = Q.front();
			Q.pop();
			for (int j = 0; j < n; ++j)
				if (A[j][u] && dis[j] == infll) {
					dis[j] = dis[u] + 1;
					Q.push(j);
				}
		}
		ans = min(ans, t + dis[n - 1]);
	}
	if (ans == infll) printf("Impossible\n");
	else printf("%lld\n", ans);
	return 0;
}

【51Nod1140】矩陣相乘結果判斷

題目連結

點選開啟連結

但是這裡需要加強一下。這裡放上加強版題目：

• 給定矩陣 \(A,B,C\) ，判斷 \(A\times B\) 是否與 \(C\) 相等
• \(T\) 組資料，\(T\le 10\)
• 每組資料，矩陣邊長 \(N(N\le 10^3)\)

題目解法

考慮向量與矩陣相乘複雜度為 \(\mathcal O(n^2)\) 。

於是就有了以下做法：隨機 5 個向量，與左邊的矩陣相乘，再與右邊矩陣相乘。判斷結果想不想等即可。

正確率極高。時間複雜度 \(\mathcal O(n^3)\) 。

程式碼

然而並沒有。

總結

• 隨機化小 trick
• 向量與矩陣相乘的複雜度的優越性

【LuoguP6569】【NOI Online #3 提高組】魔法值

題目連結

點選開啟連結

題目解法

這其實每個點的值就是由上一個狀態轉移過來的。所以很明顯可以矩陣乘法。只是這個是 \((\times,\operatorname{xor})\) 的矩陣乘法。即先對位相乘，再異或起來。

但是會發現樸素矩陣乘法是過不去的。時間複雜度為 \(\mathcal O(qn^3\log n)\) 。

考慮優化。

優化1：發現需要快速冪的矩陣是0/1矩陣，所以bitset優化即可。時間複雜度為 \(\mathcal O(\frac{qn^3\log n}{\omega})\) 。

優化2：發現需要快速冪的矩陣底數相同。所以可以預處理除矩陣的 \(2^x\) 次方。然而這沒有什麼效果。矩陣快速冪程式碼是這個樣子的（可能這份程式碼在寫的時候並不能用上，因為可能會斜掛，這是作者隨便寫的）：

int Mat(Mat x, int y) {
    Mat ret(x.n, x.m);
    ret.Init();
    for (; y; y >>= 1, x = x * x)
        if (y & 1) ret = ret * x;
    return ret;
}

你會發現預處理矩陣的 \(2^x\) 的用處只是省掉了 x=x*x 這個部分。然而 ret=ret*x 還在呢。所以得想辦法把 ret=ret*x優化一下才行。然後你作為一個智者，你發現最終需要的是 \(A\cdot B^{n}\) 這個形式的。而 \(A\) 實際是個向量。所以你想到向量乘矩陣複雜度的優越性，所以直接把 ret 換成 \(A\) 矩陣就行。時間複雜度為 \(\mathcal O(n^3\log n+qn^2 \log n)\) 。

把優化1和優化2合併起來可以得到更優的時間複雜度，\(\mathcal O(\frac{n^3\log n+qn^2 \log n}{\omega})\) 。

總結

• 向量乘矩陣複雜度的優越性。
• bitset優化0/1矩陣乘法。

程式碼

沒有。作者懶。

【CodeForces645E】Intellectual Inquiry

題目敘述

然而原題是弱化版。這裡給出加強版。

給定一個序列長度為 \(n\)，每個位置的值取 \([1,k]\)（此序列給定）現在你需要在結尾新增 \(m\) 個元素使得其本質不同的子序列儘可能多，輸出子序列數。

\(k\le 300,m\le 10^9,n\le 10^5\)

題目解法

先考慮如何求一個序列的子序列數量。

想到子序列自動機。子序列自動機是一個可以接受一個字串所有子序列的東西，一個子序列自動機上的路徑就代表一個子序列。子序列自動機上路徑數量也就是總共子序列數量。子序列自動機的一條邊上存一個數，一個結點代表一個位置。具體怎樣，就看如何構建吧！

首先一個結點即一個位置 \(i\) 向後連一條權值為 \(j\) 的邊，連到的節點是 \(i\) 後面第一個權值為 \(j\) 的位置。於是子序列自動機就是一個 DAG . 直接在上面進行 dp 就可以了。

現在貪心地考慮，想新加的一個位置 \(i\) 的貢獻是多少。設 \(f_i\) 表示位置 \(i\) 的貢獻是多少，\(l_i\) 表示上一個與 \(i\) 相同的字元的位置是哪裡。那麼有 \(f_i=\sum_{j=l_i+1}^{i-1}f_j\) 。所以新加的位置應該滿足上一次出現位置儘量靠前，這樣能讓當前位置儘量大。

然後試幾個數就可以發現後面填的數是有周期規律的。所以考慮用 \(f_i\) 的遞推式快速推出 \(f\) 數列的第 \(n\) 項。

總結

• 子序列自動機求子序列數量

程式碼

並沒有。

【LuoguT134222】【YsOI團隊主の神仙題賽】神仙題

這題我也不知道哪裡來的，但是題還不錯，所以講講。

注：這道題中 \(\binom{x}{y}=\mathrm{C}_{x}^{y}\) ，可能有些人不知道這個記號。

題目連結

點選開啟連結

題目解法

推式子。對於一個固定的 \(m\) 我們設 \(f_i\) 表示 \(n=i\) 時原式的值。

\[\begin{aligned} f_n&=\sum_{i\le n-i\cdot m} \binom{n-i\cdot m}{i}\\ &=\sum_{x+my=n}\binom{x}{y}\\ &=\sum_{x+my=n}\left(\binom{x-1}{y-1} + \binom{x-1}{y}\right)\\ &=\sum_{x+my=n-m-1}\binom{x}{y}+\sum_{x+my=n-1}\binom{x}{y}\\ &=f_{n-m-1}+f_{n-1} \end{aligned} \]

所以矩陣快速冪遞推轉移即可。

總結

• 使用 \(\binom{x}{y}=\binom{x-1}{y-1}+\binom{x-1}{y}\) 把組合數求和轉化為遞推式。

程式碼

莫得程式碼。

【LOJ2143】「SHOI2017」組合數問題

題目連結

點選開啟連結

題解

待補充。

總結

待補充。

程式碼

莫得程式碼。