Luogu P5279 [ZJOI2019]麻将

ZJOI2019神题，间接送我退役的神题233

考场上由于T2写挂去写爆搜的时候已经没多少时间了，所以就写挂了233

这里不多废话直接开始讲正解吧，我们把算法分成两部分

---

1.建一个“胡牌自动机”

首先我们发现这题不能转化为一般DP问题求解的最大瓶颈就是因为它的状态很诡异

但是我们细细一想，形如({1,1,1,2,3,4,5,6,7,8,9,10,12,12})和({3,3,3,5,6,7,9,10,11,14,15,16,20,20})的本质其实是一样的（都是三个顺子+一个刻子+一个雀头的形式，且相对大小分布相同）

那么考虑到一共只有(13)张牌，那么最后的有效集合数必然很少，如果我们可以构造出这些集合那么可以帮助DP了

如果你知道DP套DP的姿势，那么就会套路地明白内层的DP应该用自动机的转移形式，所以我们大致的目标就出现了：

建出一个自动机，其中每一个节点都是一个状态，然后每个节点向子节点连的边就是摸牌之后能转移到的状态

所以我们就先考虑怎么建这个自动机，还是考虑DP

我们设(f_{0/1,i,j,k})表示目前处理完了前(i)种牌，还剩下(j)组((i-1,i))以及(k)张(i)，且是否（用(0/1)表示）存在雀头时最多的面子数

由于当(j,kge 3)时可以直接用(i-1,i)组成刻子，因此我们状态中的(j,kle 2)

那么我们就开一个(3 imes3)的矩阵表示状态，单个转移的时候枚举用于拼面子，与((i-1,i))拼顺子，以及用来保留的个数，然后将多余的拿来拼刻子即可

然后考虑每次多出一个数值，我们直接讨论是否要留出一个雀头，分情况转移即可

我们发现这样并没有考虑七对子的情况，这个没关系，我们直接记录雀头的个数，特判了即可

因此剩下终止状态的判断就很简单了，直接把能胡的点作为终止节点结束即可

具体构造的过程可以用一个map来去重，然后用BFS来扩展状态，这个具体看代码

---

2.期望DP

先说一句，前面由于自动机对于任意数据构造相同，所以总点数是固定的(2092)，如前言所述不大

那么考虑DP求解最后的问题，我们用一个经典套路，将求(ge)的化为求(>)然后最后加上等于的情况即可

具体到这道题上，其实就是求出(i)轮后不胡的方案数然后乘上贡献，最后的总答案加(1)即可

那么大致的DP方程就有了，我们设(f_{i,j,k})表示选了前(i)种牌，一共用了(j)张，此时位于自动机上(k)号点的方案数

那么转移其实很简单，枚举第(i)种牌选的张数(t)，(f_{i,j,k})向(f_{i+1,j+t,node_k.son_{a_i+t}})转移，注意不要忘记乘上组合数(C_{4-a_i}^t)

求出(f)数组后最后的答案统计就十分简单了，令(m=4n-13)：

[ans=sum_{i=1}^m frac{f_icdot i!cdot(m-i)!}{m!}+1 ]

这个很好理解，因为前后都可以随便选，因此不再赘述

综上，我们得到了一个复杂度上界为(2092cdot n^2)的优秀做法，因为常数很小且跑不满因此足以通过

CODE

#include<cstdio>
#include<cstring>
#include<map>
#define RI register int
#define CI const int&
using namespace std;
const int N=105,MX=2100,mod=998244353;
inline void maxer(int& x,CI y)
{
    if (y>x) x=y;
}
inline int min(CI a,CI b)
{
    return a<b?a:b;
}
struct Matrix
{
    int mat[3][3];
    inline int* operator [] (CI x) { return mat[x]; }
    inline Matrix(void)
    {
        memset(mat,-1,sizeof(mat));
    }
    friend inline bool operator != (Matrix A,Matrix B)
    {
        for (RI i=0;i<3;++i) for (RI j=0;j<3;++j)
        if (A[i][j]!=B[i][j]) return 1; return 0;
    }
    friend inline bool operator < (Matrix A,Matrix B)
    {
        for (RI i=0;i<3;++i) for (RI j=0;j<3;++j)
        if (A[i][j]!=B[i][j]) return A[i][j]<B[i][j];
    }
    inline bool Com_Hu(void)
    {
        for (RI i=0;i<3;++i) for (RI j=0;j<3;++j)
        if (mat[i][j]>=4) return 1; return 0;
    }
    inline void flush(Matrix pre,CI num)
    {
        for (RI i=0;i<3;++i) for (RI j=0;j<3;++j)
        if (~pre[i][j]) for (RI k=0;k<3&&i+j+k<=num;++k)
        maxer(mat[j][k],min(i+pre[i][j]+(num-i-j-k)/3,4));
    }
};
struct Hu_Auto_Node
{
    Matrix p[2]; int cur,ch[5];
    inline Hu_Auto_Node(void)
    {
        memset(ch,0,sizeof(ch)); cur=0; p[0]=p[1]=Matrix();
    }
    inline bool is_Hu(void)
    {
        if (cur>=7) return 1; return p[1].Com_Hu();
    }
    inline void Hu(void)
    {
        memset(ch,0,sizeof(ch)); cur=-1; p[0]=p[1]=Matrix();
    }
    friend inline bool operator < (Hu_Auto_Node A,Hu_Auto_Node B)
    {
        if (A.cur!=B.cur) return A.cur<B.cur;
        if (A.p[0]!=B.p[0]) return A.p[0]<B.p[0];
        if (A.p[1]!=B.p[1]) return A.p[1]<B.p[1]; return 0;
    }
    friend inline Hu_Auto_Node operator + (Hu_Auto_Node A,CI num)
    {
        if (A.is_Hu()) return A.Hu(),A; Hu_Auto_Node s;
        s.p[0].flush(A.p[0],num); s.p[1].flush(A.p[1],num);
        if (num>=2) s.p[1].flush(A.p[0],num-2);
        s.cur=A.cur+(num>=2); if (s.is_Hu()) s.Hu(); return s;
    }
};
class Hu_Automation
{
    private:
        map <Hu_Auto_Node,int> Hash;
        inline void expand(CI id)
        {
            for (RI i=0;i<=4;++i)
            {
                Hu_Auto_Node son=node[id]+i;
                if (!Hash.count(son)) node[Hash[son]=++tot]=son;
                node[id].ch[i]=Hash[son];
            }
        }
    public:
        Hu_Auto_Node node[2100]; int tot;
        inline Hu_Automation(void)
        {
            node[1].p[0][0][0]=0; node[tot=2].cur=-1;
            Hash[node[1]]=1; Hash[node[2]]=2; expand(1);
            for (RI i=3;i<=tot;++i) expand(i);
        }
        /*inline void check(void)
        {
            printf("%d
",tot); for (RI i=1;i<=tot;++i,putchar('
'))
            for (RI j=0;j<=4;++j) printf("%d ",node[i].ch[j]);
        }*/
}HA;
int f[2][N<<2][MX],a[N],fact[N<<2],inv[N<<2],n,m,x,y,nw,ans;
inline void inc(int& x,CI y)
{
    if ((x+=y)>=mod) x-=mod;
}
inline int quick_pow(int x,int p=mod-2,int mul=1)
{
    for (;p;p>>=1,x=1LL*x*x%mod) if (p&1) mul=1LL*mul*x%mod; return mul;
}
inline void init(CI n)
{
    RI i; for (fact[0]=i=1;i<=n;++i) fact[i]=1LL*fact[i-1]*i%mod;
    for (inv[n]=quick_pow(fact[n]),i=n-1;~i;--i) inv[i]=1LL*inv[i+1]*(i+1)%mod;
}
inline int C(CI n,CI m)
{
    return 1LL*fact[n]*inv[m]%mod*inv[n-m]%mod;
}
inline void calc(CI x,CI y)
{
    inc(ans,1LL*f[nw][x][y]*fact[x]%mod*fact[m-x]%mod);
}
int main()
{
    RI i,j,k,t; for (scanf("%d",&n),i=1;i<=13;++i) scanf("%d%d",&x,&y),++a[x];
    for (init(m=(n<<2)-13),f[0][0][1]=i=1;i<=n;++i)
    for (nw=i&1,memset(f[nw],0,sizeof(f[nw])),j=m;~j;--j)
    for (k=1;k<=HA.tot;++k) if (f[nw^1][j][k]) for (t=0;t<=4-a[i];++t)
    inc(f[nw][j+t][HA.node[k].ch[a[i]+t]],1LL*f[nw^1][j][k]*C(4-a[i],t)%mod);
    for (nw=n&1,i=1;i<=m;++i) for (calc(i,1),j=3;j<=HA.tot;++j) calc(i,j);
    return printf("%d",(1LL*ans*inv[m]%mod+1)%mod),0;
}