Problem
Bessie 最近参加了一场 USACO 竞赛，遇到了以下问题。当然 Bessie 知道怎么做。那你呢？
考虑一个仅由范围在 \(1 \ldots K\)（\(1 \leq K \leq 20\)）之间的整数组成的长为 \(N\) 的序列 \(A_1,A_2, \ldots ,A_N\)（\(1 \leq N \leq 5 \times 10^4\)）。给定 \(Q\)（ \(1 \leq Q \leq 2 \times 10^5\) ）个形式为 \([L_i,R_i]\)（\(1 \leq L_i \leq R_i \leq N\)）的询问。对于每个询问，计算 \(A_{L_i},A_{L_i+1}, \ldots ,A_{R_i}\) 中不下降子序列的数量模 \(10^9+7\) 的余数。
\(A_L,\ldots ,A_R\) 的一个不下降子序列是一组索引 （\(j_1,j_2, \ldots ,j_x\)），满足 \(L\le j_1<j_2<\ldots<j_x\le R\) 以及 \(A_{j_1}\le A_{j_2}\le \ldots \le A_{j_x}\)。确保你考虑了空子序列！
Sol
\(f_{i,j}\) 为考虑前 \(i\) 个，子序列结尾值为 \(j\) 方案数。有转移 \(f_{i, a_i} = f_{i-1, a_i} + \sum_{j \le a_i} f_{i-1, j}\)，对于 \(j \ne a_i\) 有 \(f_{i, j} = f_{i - 1, j}\)。\(j\) 这一维很小，考虑矩阵优化。转移矩阵 \(A\) 等于：
推一下可以简单的得到逆矩阵。即单位矩阵第 \(a_i\) 行前面 \(a_i-1\) 个都是 \(-\frac12\) 然后剩下的是单位矩阵。直接做矩阵的前缀（逆）积是 \(nk^3\) 无法接受。发现有值的位置仅 \(O(k)\) 个，可以优化到 \(O(nk^2)\)。查询可以用 base 来乘是 \(k^2\)，可以前缀和做到 \(O(k)\)。
假如没有逆矩阵怎么做？
线段树维护矩阵是 \(O(nk^3)-O(k^2\log n)\)。猫树是 \(O(nk^3\log)-O(qk^2)\)。若使用分块减少预处理量，散点转移前缀和，也难以通过。
换个思路，直接使用线段树维护。节点存储开始为 \(x\)，结束为 \(y\) 的答案。使用前缀和合并减到 \(O(k^2)\)。然后再使用猫树。空间 \(O(nk^2)\) 可以通过离线逐块处理（即变成分治）做到 \(O(k^2\log n)\)。时间复杂度 \(O(nk^2\log n + qk^2)\)。
Code
代码是写的逆矩阵的。远古代码，还请见谅。
#include<bits/stdc++.h>
#define pb push_back
#define ios ios::sync_with_stdio(0); cin.tie(0);
using namespace std;
typedef long long ll;
typedef pair<int, int> pii;
const int inf = 1e9 + 7, infi32 = 2147483647, mod = 1e9 + 7;
const ll INF = 1e18;const int N = 5e4 + 10, MAXK = 20 + 1, inv2 = (mod + 1) / 2;
int n, K, Q;
int a[N], sp[N][MAXK];
struct Matrix {int a[MAXK][MAXK];
} pre[N], ipre[N];int main() {ioscin >> n >> K;for (int i = 1; i <= n; i++) cin >> a[i];for (int i = 0; i <= K; i++)pre[0].a[i][i] = ipre[0].a[i][i] = 1;for (int i = 1; i <= n; i++) {pre[i] = pre[i - 1]; ipre[i] = ipre[i - 1];for (int j = 0; j <= K; j++)for (int k = 0; k <= a[i]; k++)pre[i].a[j][a[i]] = (pre[i].a[j][a[i]] + pre[i - 1].a[j][k]) % mod;for (int j = 0; j <= a[i]; j++)for (int k = 0; k <= K; k++)ipre[i].a[j][k] = (ipre[i].a[j][k] + 1ll * ipre[i - 1].a[a[i]][k] * (mod - inv2) % mod) % mod;for (int j = 0; j <= K; j++)for (int k = 0; k <= K; k++)sp[i][j] = (sp[i][j] + pre[i].a[j][k]) % mod;}cin >> Q;while (Q--) {int l, r, ans = 0; cin >> l >> r;for (int i = 0; i <= K; i++) ans = (ans + 1ll * ipre[l - 1].a[0][i] * sp[r][i] % mod) % mod;cout << ans << "\n";}return 0;
}