AC 自动机 DP 模板
以及矩阵快速幂优化线性递推 DP。
本篇主要讲矩阵快速幂优化线性递推 DP。
题目
NOIP3-T2 佩绣-ACGO
题意
给 \(m\) 个由小写英文字母构成的字符串 \(S_i\),确定在长度不超过 \(n\) 的字符串 \(T\) 满足 \(T\) 的任意子串都不等于 \(S_i\) 的数量。
答案对 \(998244353\) 取模。
输入
op n m
S_1
...
S_m
输出
ans
数据范围
- \(1 \le op \le 10\)
- \(1 \le n \le 10^{18}\)
- \(1 \le m \le 200\)
- \(\displaystyle \sum_{i=1}^{m} |S_i| \le 200\)
Subtask 1:\(O(n \cdot sz)\)
记 trie 树上状态数为 \(sz\)。
由 AC 自动机可以求出哪些状态是合法的。
考虑 DP,定义 \(dp_{len}[v]\) 表示当长度为 \(len\) 在 AC 自动机匹配状态为 \(v\) 的方案数。
\(v\) 的最后一个字母为 \(c\):
\[dp_{\,len}[v] = \sum_{u \notin \text{bad}} dp_{\,len-1}[u], \quad v = \text{trie}[u][c],\ c \in [0, 26)
\]
\[ans = \sum_{len=1}^{n} \sum_{u \notin \text{bad}} dp_{\,len}[u]
\]
注意:当 AC 自动机求失配指针时遇到不存在的状态 \(trie_{u,c}\) 时,请将 \(trie_{u,c} \gets trie_{fail_u, c}\)。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=998244353;
const ll maxn=1e6+5;
int trie[maxn][30],fail[maxn];
int cnt;
int bad[maxn];
void insert(string s) {ll now=0;for (char c:s) {ll ch=c-'a';if (!trie[now][ch]) trie[now][ch]=++cnt;now=trie[now][ch];}bad[now]=1;
}
void build() {queue<ll> q;for (ll i=0;i<26;i++) if (trie[0][i]) q.push(trie[0][i]);while (!q.empty()) {ll now=q.front();q.pop();bad[now]|=bad[fail[now]];for (ll i=0;i<26;i++) {if (trie[now][i]) {fail[trie[now][i]]=trie[fail[now]][i];q.push(trie[now][i]);}else trie[now][i]=trie[fail[now]][i];}}
}
void solve() {ll op,n,m;cin>>op>>n>>m;for (ll i=1;i<=m;i++) {string s;cin>>s;insert(s);}build();ll sz=cnt;vector<ll> dp(sz+1,0);dp[0]=1;ll ans=0;for (ll len=1;len<=n;len++) {vector<ll> ndp(sz+1,0);for (ll u=0;u<=sz;u++) {if (!dp[u]||bad[u]) continue;for (ll c=0;c<26;c++) {ll v=trie[u][c];if (!bad[v]) ndp[v]=(ndp[v]+dp[u])%mod;}}dp=ndp;for (ll u=0;u<=sz;u++) ans=(ans+dp[u])%mod;}cout<<ans<<"\n";
}
int main() {ios::sync_with_stdio(0);ll t=1;while (t--) solve();
}
Subtask 2:矩阵快速幂优化
\[dp_{\,len}[v] = \sum_{u \notin \text{bad}} dp_{\,len-1}[u], \quad v = \text{trie}[u][c],\ c \in [0, 26)
\]
\[ans = \sum_{len=1}^{n} \sum_{u \notin \text{bad}} dp_{\,len}[u]
\]
每次转移只依赖上一轮,是线性递推。
转移矩阵 \(T \in \mathbb{Z}^{sz \times sz}\):
\[T[v][u] = \begin{cases}
1 & \text{从 } u \text{ 走一步到 } v \text{ 合法} \\
0 & \text{否则}
\end{cases}
\]
其实 \(T[v][u]\) 就是 trie 树上从 \(u\) 到 \(v\) 的边数(\(0\) 或 \(1\)):
\[dp_{\,len}[v] = \sum_{u=0}^{sz-1} T[v][u] \cdot dp_{\,len-1}[u]
\]
具体矩阵转移:
\[\begin{bmatrix}
dp_{\,len}[0] \\
dp_{\,len}[1] \\
\vdots \\
dp_{\,len}[sz-1]
\end{bmatrix}
=
\begin{bmatrix}
T[0][0] & T[0][1] & \cdots & T[0][sz-1] \\
T[1][0] & T[1][1] & \cdots & T[1][sz-1] \\
\vdots & \vdots & \ddots & \vdots \\
T[sz-1][0] & T[sz-1][1] & \cdots & T[sz-1][sz-1]
\end{bmatrix}
\cdot
\begin{bmatrix}
dp_{\,len-1}[0] \\
dp_{\,len-1}[1] \\
\vdots \\
dp_{\,len-1}[sz-1]
\end{bmatrix}
\]
设状态向量 \(V_{len}\) 长度为 \(sz\):
\[V_{len}[u] = dp_{\,len}[u]
\]
于是:
\[V_{len} = T \cdot V_{len-1}
\]
\[V_{len} = T^{\,len} \cdot V_0
\]
暴力 DP 每轮累加求出答案,而矩阵版需要一次算出。
扩充矩阵,多一行一列。
设扩充矩阵 \(M \in \mathbb{Z}^{(sz+1) \times (sz+1)}\):
\[M[i][j] =
\begin{cases}
T[i][j] & (0 \le i,j < sz) \\[4pt]
\displaystyle\sum_{v} T[v][j] & (i = sz,\ 0 \le j < sz) \\[4pt]
1 & (i = sz,\ j = sz) \\[4pt]
0 & \text{否则}
\end{cases}
\]
- \(M[0..sz-1][0..sz-1] = T\):保持原转移。
- \(M[sz][j] = \sum_v T[v][j]\):\(j\) 的所有出度,即当前步 \(j\) 产生的合法串数,累加到
ans行。 - \(M[sz][sz] = 1\):保留之前累加的
ans。
具体转移:
\[\begin{bmatrix}
dp_{\,len}[0] \\
dp_{\,len}[1] \\
\vdots \\
dp_{\,len}[sz-1] \\
\hline
ans_{\,len}
\end{bmatrix}
=
\begin{bmatrix}
T[0][0] & \cdots & T[0][sz-1] & 0 \\
\vdots & \ddots & \vdots & \vdots \\
T[sz-1][0] & \cdots & T[sz-1][sz-1] & 0 \\
\hline
\sum_v T[v][0] & \cdots & \sum_v T[v][sz-1] & 1
\end{bmatrix}
\cdot
\begin{bmatrix}
dp_{\,len-1}[0] \\
\vdots \\
dp_{\,len-1}[sz-1] \\
\hline
ans_{\,len-1}
\end{bmatrix}
\]
最后一行展开:
\[ans_{\,len} = \underbrace{\sum_{u=0}^{sz-1} \left( \sum_{v} T[v][u] \right) \cdot dp_{\,len-1}[u]}_{\text{新贡献 } \sum dp_{\,len}} + \underbrace{ans_{\,len-1}}_{\text{旧和}}
\]
初始值:
\[V'_0 =
\begin{bmatrix}
1 \\
0 \\
\vdots \\
0 \\
\hline
0
\end{bmatrix}
\quad
\text{(空串在根节点 $0$,`ans` 初始为 $0$)}
\]
\(M^n \cdot V'_0\) 的最后一行的值就是答案。
最后矩阵快速幂优化即可从 \(O(n \cdot sz)\) 变为 \(O(sz^3 \log n)\)。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=998244353;
const ll maxn=1e6+5;
int trie[maxn][30],fail[maxn];
int cnt;
int bad[maxn];
void insert(string s) {ll now=0;for (char c:s) {ll ch=c-'a';if (!trie[now][ch]) trie[now][ch]=++cnt;now=trie[now][ch];}bad[now]=1;
}
void build() {queue<ll> q;for (ll i=0;i<26;i++) if (trie[0][i]) q.push(trie[0][i]);while (!q.empty()) {ll now=q.front();q.pop();bad[now]|=bad[fail[now]];for (ll i=0;i<26;i++) {if (trie[now][i]) {fail[trie[now][i]]=trie[fail[now]][i];q.push(trie[now][i]);}else trie[now][i]=trie[fail[now]][i];}}
}
struct nd {ll n;vector<vector<ll>> a;nd(ll _n) {n=_n;a.resize(n,vector<ll>(n,0));}nd operator * (const nd &b) const {nd res(n);for (ll i=0;i<n;i++)for (ll k=0;k<n;k++) if (a[i][k])for (ll j=0;j<n;j++) if (b.a[k][j])res.a[i][j]=(res.a[i][j]+a[i][k]*b.a[k][j])%mod;return res;}
};
nd qpow(nd a,ll b) {nd res(a.n);for (ll i=0;i<a.n;i++) res.a[i][i]=1;while (b) {if (b&1) res=res*a;a=a*a;b>>=1;}return res;
}
void solve() {ll op,n,m;cin>>op>>n>>m;for (ll i=1;i<=m;i++) {string s;cin>>s;insert(s);}build();ll sz=cnt+1;nd M(sz+1);for (ll u=0;u<sz;u++) {if (bad[u]) continue;ll sum=0;for (ll c=0;c<26;c++) {ll v=trie[u][c];if (!bad[v]) {M.a[v][u]++;sum++;}}M.a[sz][u]=sum;}M.a[sz][sz]=1;M=qpow(M,n);cout<<M.a[sz][0]<<endl;
}
int main() {ios::sync_with_stdio(0);ll t=1;while (t--) solve();
}
