昨晚搓了一局雀魂，牌运还不错；今天和同学聊天的时候提到了，然后顺便翻出了这道题和隔壁 GXOI/GZOI 那道麻将题。
然后顺便研究了一下这道题，其实也不是特别难写，不过想出来应该也很难。
首先我们考虑如何判定一堆牌的子集能否和牌。
注意到我们并不关心牌之间的顺序，故可以将这堆牌表示为「第 \(i\) 种牌有 \(a_i\) 种」的形式。
不难想到 DP，设 \(f_{i,0/1,j,k}\) 表示考虑到第 \(i\) 种牌，是否预留过雀头，预留了 \(j\) 组 \((i-1,i)\) 和 \(k\) 张 \(i\) 的情况下能凑出的最多的面子数。
并且因为 \(j \ge 3\) 或 \(k \ge 3\) 都可以算为面子，所以 \(j < 3,k < 3\)。
如果我们要多次判定但只能做一次 DP 呢？
可以考虑将这个 DP 改成自动机。具体地，考虑把 \(f\) 数组删去第一维得到的数组作为一个自动机上的结点，DP 转移对应自动机上的转移。
在构造的时候，不考虑 \(i\) 的影响，因为最后无论如何都一定会得到和牌的状态。而得到和牌状态时则不能再改变状态，否则便会无限进行下去。
这样可以把相同结点合并，经实测，最多会有 \(2092\) 个结点。
判定的时候，只需要在自动机上跑边即可。
那么此题显然可以考虑在自动机上再设计一个 DP，达到 DP 套 DP 的效果。
考虑设 \(g_{i,j,k}\) 表示考虑到第 \(i\) 种牌，已经摸了 \(j\) 张牌，走到自动机上的 \(k\) 号结点有多少种方案。
可以容易地从这个 DP 得到 \(g'_i\) 表示摸 \(i\) 张牌都不和牌的方案数，由此便可知答案为 \[
1+\frac1{(4n-13)!}\sum\limits_{i=1}^{4n-13} g'_i i!(4n-13-i)!
\]
因为如果某种排列下摸 \(i+1\) 张牌才能和牌而 \(i\) 张牌则不行，那么对 \(i-1,i-2,\ldots\) 也会有贡献，总共就是刚好 \(i\) 次贡献。 并且确定了排列的前 \(i\) 张牌是哪些牌，顺序无所谓，后面的牌顺序也无所谓，所以要乘阶乘。
而题目问的是第一张能和牌的牌，所以要加一。
所以好像还没说 \(g\) 的转移？
也蛮简单的吧，注意枚举取多少第 \(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
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
using namespace std;
const int N = 100;
const int M = 2092;
const int mod = 998244353;
inline int fpow(int a,int b)
{
int ret = 1;
for(;b;b >>= 1)
(b & 1) && (ret = (long long)ret * a % mod),a = (long long)a * a % mod;
return ret;
}
int n,a[N + 5];
struct Matrix
{
int a[3][3];
inline int *operator[](const int &x)
{
return a[x];
}
inline const int *operator[](const int &x) const
{
return a[x];
}
inline Matrix()
{
for(register int i = 0;i < 3;++i)
for(register int j = 0;j < 3;++j)
a[i][j] = -1;
}
inline bool operator<(const Matrix &o) const
{
for(register int i = 0;i < 3;++i)
for(register int j = 0;j < 3;++j)
if(a[i][j] ^ o[i][j])
return a[i][j] < o[i][j];
return 0;
}
inline bool operator!=(const Matrix &o) const
{
for(register int i = 0;i < 3;++i)
for(register int j = 0;j < 3;++j)
if(a[i][j] ^ o[i][j])
return 1;
return 0;
}
inline void trans(Matrix &o,int x) const
{
for(register int i = 0;i < 3;++i)
for(register int j = 0;j < 3;++j)
if(~a[i][j])
for(register int k = 0;k < 3 && i + j + k <= x;++k)
o[j][k] = max(o[j][k],a[i][j] + i + (x - i - j - k) / 3),
o[j][k] = min(o[j][k],4);
}
inline bool hou() const
{
for(register int i = 0;i < 3;++i)
for(register int j = 0;j < 3;++j)
if(a[i][j] == 4)
return 1;
return 0;
}
};
struct node
{
Matrix jantou[2];
int toitsu;
inline bool operator<(const node &o) const
{
if(toitsu ^ o.toitsu)
return toitsu < o.toitsu;
if(jantou[0] != o.jantou[0])
return jantou[0] < o.jantou[0];
if(jantou[1] != o.jantou[1])
return jantou[1] < o.jantou[1];
return 0;
}
inline void clear()
{
jantou[0] = jantou[1] = Matrix(),toitsu = 0;
}
inline void begin()
{
clear(),jantou[0][0][0] = 0;
}
inline void end()
{
clear(),toitsu = -1;
}
inline node()
{
clear();
}
inline node(int op)
{
op ? end() : begin();
}
inline bool hou() const
{
if(toitsu >= 7)
return 1;
if(jantou[1].hou())
return 1;
return 0;
}
inline node operator+(int x) const
{
if(toitsu == -1)
return *this;
node ret;
x >= 2 && (jantou[0].trans(ret.jantou[1],x - 2),1);
jantou[0].trans(ret.jantou[0],x);
jantou[1].trans(ret.jantou[1],x);
ret.toitsu = toitsu + (x >= 2);
ret.hou() && (ret.end(),1);
return ret;
}
} q[M + 5];
int head,tail;
int tot,e[M + 5][5],ed;
map<node,int> id;
inline void bfs()
{
id[q[++tail] = node(0)] = ++tot;
for(register node p,dir;head < tail;)
{
p = q[++head];
int x = id[p];
for(register int i = 0;i <= 4;++i)
{
dir = p + i;
if(id.count(dir))
e[x][i] = id[dir];
else
e[x][i] = id[dir] = ++tot,q[++tail] = dir;
}
}
ed = id[node(1)];
}
int fac[(N << 2) + 5],ifac[(N << 2) + 5];
int f[2][(N << 2) + 5][M + 5],ans;
int main()
{
bfs();
scanf("%d",&n),fac[0] = 1;
for(register int i = 1;i <= (n << 2);++i)
fac[i] = (long long)fac[i - 1] * i % mod;
ifac[n << 2] = fpow(fac[n << 2],mod - 2);
for(register int i = n << 2;i;--i)
ifac[i - 1] = (long long)ifac[i] * i % mod;
int x;
for(register int i = 1;i <= 13;++i)
scanf("%d%*d",&x),++a[x];
f[0][0][id[node(0)]] = 1;
for(register int i = 1;i <= n;++i)
{
memset(f[i & 1],0,sizeof f[i & 1]);
for(register int j = 0;j <= (n << 2) - 13;++j)
for(register int k = 1;k <= tot;++k)
if(f[(i & 1) ^ 1][j][k])
for(register int l = 0;l <= 4 - a[i] && i + l <= (n << 2) - 13;++l)
f[i & 1][j + l][e[k][a[i] + l]] = (f[i & 1][j + l][e[k][a[i] + l]] + (long long)f[(i & 1) ^ 1][j][k] * fac[4 - a[i]] % mod * ifac[l] % mod * ifac[4 - a[i] - l] % mod) % mod;
}
for(register int i = 1;i <= (n << 2) - 13;++i)
for(register int j = 1;j <= tot;++j)
(j ^ ed) && (ans = (ans + (long long)f[n & 1][i][j] * fac[i] % mod * fac[(n << 2) - 13 - i] % mod) % mod);
ans = ((long long)ans * ifac[(n << 2) - 13] % mod + 1) % mod;
printf("%d\n",ans);
}