首先看到这个数据范围限制，我们不难猜到是状压 \(dp\)，首先就猜测复杂度应该是 \(O(nmS)\)，\(S\) 是状压大小，和 \(m\) 有关的某个非多项式的级数。
错误的思考过程
然后我们就很快能找到一个状压 \(dp\) 的方法。我们可以记录当前一列的黑白染色情况和连通块的集合划分。粗略估计是 \(\sum_{i=0}^{m}\binom{m}{i}Bell_iBell_{m-i}\) 的，很寄就是了。
接下来我们考虑如何节省无用的状态。首先，我们发现，在当前列相邻的两个同色的格子是一定是相同连通块的。那么我们可以把状态中相邻的相同颜色合并。因此，每个状态中，最多只有四个相同颜色的部分需要划分。而 \(Bell_4\) 是很小的，只有 \(15\)。而且取到 \(4\) 的条件是非常苛刻的。如果我们仔细计算一下，我们会发现合法的状态一共就只有 \(3060\) 种。
由此我们就欣喜若狂了，我们的复杂度似乎可以做到 \(3.06\cdot 10^7\) 次运算，这是非常优秀的。但是一旦我们冷静下来就会发现，这个东西的转移非常困难。
首先，这些序列和状态之间没有好的哈希关系，这就需要我们在程序的一开始暴力枚举从而对状态进行编码。其次，我们为了 \(Bell_4\) 分开了白和黑，这就导致将其重构成序列再压缩回集合划分的步骤繁琐。
其次，我们每次都需要重新构造序列，然后再把集合划分序列转成哈希值存进状态里。这个过程是 \(O(m)\) 的。并且转移很多步骤都非常苛刻，难以统一，需要分类讨论很多。
总而言之，一开始不顾转移对状态的过度精简，导致这个程序怀着美好的梦想但几乎无法实现。下面贴一下按照这个思路实现出来的代码，\(Bug\) 连篇，可能什么时候码力上升了会试着回来调一下。但现在还是算了。
ll n,m,P;
int bel[5]={1,1,2,5,15};
int stac;
int divc;
int divt[4][4][4][4];
int stat[256][15][15];
array<int,4> divh[15];
array<int,3> stah[3060];
int cb[2][256];
namespace init{
inline void init_div(){
array<int,4>ori;
vt<array<int,4> >res;
function<void(int,int,array<int,4>)>dfs=[&](int x,int y,array<int,4>cur) {
if(x==4)res.pb(cur);
else rep(j,0,y+1)cur[x]=j,dfs(x+1,max(j,y),cur);
};
dfs(0,-1,ori);
sort(res.begin(),res.end(),[](array<int,4>a,array<int,4>b){
if(a[3]!=b[3])return a[3]<b[3];
if(a[2]!=b[2])return a[2]<b[2];
if(a[1]!=b[1])return a[1]<b[1];
return a[0]<b[0];
});
for(auto a:res){
divt[a[0]][a[1]][a[2]][a[3]]=divc;
divh[divc++]=a;
}
}
inline void init_sta(){
for(int msk=0;msk<(1<<m);msk++){
cb[msk&1][msk]++;
rep(j,1,m-1)cb[msk>>j&1][msk]+=((msk>>j&1)!=(msk>>(j-1)&1));
for(int i=0;i<bel[cb[0][msk]];i++){
for(int j=0;j<bel[cb[1][msk]];j++){
stat[msk][i][j]=stac;
stah[stac++]={msk,i,j};
}
}
}
}
};
using namespace init;
int dp[3060][4],f[3060][4],color[8],cur[8],num[8];
int cnt[2][5];
inline void restruct(int &a,int &b){
array<int,4>A={0,0,0,0},B={0,0,0,0};
int ca=0,cb=0,cA=0,cB=0;
rd(i,2)rd(j,5)cnt[i][j]=0;
rep(i,0,7)if(i==0||color[i]!=color[i-1]){
if(color[i]){
if(!cnt[1][cur[i]])cnt[1][cur[i]]=ca++;
A[cA++]=cnt[1][cur[i]];
}else{
if(!cnt[0][cur[i]])cnt[0][cur[i]]=cb++;
B[cB++]=cnt[0][cur[i]];
}
}
a=divt[A[0]][A[1]][A[2]][A[3]];
b=divt[B[0]][B[1]][B[2]][B[3]];
}
inline void add(int &x,int y){
x=(x+y)%P;
}
inline void solve(int y){
//现在要填第y列
rd(msk,(1<<m))rd(bl,bel[cb[0][msk]]-1)rd(wh,bel[cb[1][msk]]-1)rd(sa,4){
if(!f[stat[msk][bl][wh]][sa])continue;
int res=f[stat[msk][bl][wh]][sa];
//现在安排的是这样的。
//分类现在是0还是1
int n[2]={0,0},cnt=0;
rd(i,m){
if(msk>>i&1){
color[i]=1,num[i]=n[1];
cur[i]=divh[bl][n[1]];
}else{
color[i]=0,num[i]=n[0];
cur[i]=divh[wh][n[0]];
}
if(i!=7&&(msk>>i&1)!=(msk>>(i+1)&1))n[msk>>i&1]++;
}
rd(i,m)if(color[i]==color[y]&&cur[i]==cur[y])cnt++;
if(y!=0&&sa!=cur[y-1])cnt++;
//左边是什么
if(color[y]){
if(y==0){
add(dp[stat[msk][bl][wh]][sa],res);
if(y!=8&&color[y+1]){
add(dp[stat[msk][bl][wh]][sa],res);
}else if(cnt!=1){
int nb,nw;
color[y]=0,cur[y]=4;
restruct(nb,nw);
add(dp[stat[msk^(1<<y)][nb][nw]][(y==8?cur[0]:color[y+1])],res);
}
}else if(y==0||!color[y-1]){
add(dp[stat[msk][bl][wh]][sa],res);
if(y!=8&&color[y+1]){
add(dp[stat[msk][bl][wh]][sa],res);
}else if(cnt!=1){
int nb,nw;
color[y]=0,cur[y]=cur[y-1];
restruct(nb,nw);
add(dp[stat[msk^(1<<y)][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}
}else{//上面是1
if(sa==cur[y-1]){//上面连通块一样
//填1
add(dp[stat[msk][bl][wh]][sa],res);
//填0
int nb,nw;
color[y]=0,cur[y]=4;
restruct(nb,nw);
add(dp[stat[msk^(1<<y)][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}else{//上面连通块不一样
//填0
if(cnt!=1){
int nb,nw;
color[y]=0,cur[y]=4;
restruct(nb,nw);
add(dp[stat[msk^(1<<y)][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}
//填1
int nb,nw,dc=cur[y-1],df=cur[y];
rd(j,8)if(color[j]&&cur[j]==df)cur[j]=dc;
restruct(nb,nw);
add(dp[stat[msk][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}
}
}else{//左边是0
if(y==0){
add(dp[stat[msk][bl][wh]][sa],res);
if(y!=8&&!color[y+1]){
add(dp[stat[msk][bl][wh]][sa],res);
}else if(cnt!=1){
int nb,nw;
color[y]=1,cur[y]=4;
restruct(nb,nw);
add(dp[stat[msk^(1<<y)][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}
}else if(color[y-1]){
add(dp[stat[msk][bl][wh]][sa],res);
if(y!=8&&!color[y+1]){
add(dp[stat[msk][bl][wh]][sa],res);
}else if(cnt!=1){
int nb,nw;
color[y]=1,cur[y]=cur[y-1];
restruct(nb,nw);
add(dp[stat[msk^(1<<y)][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}
}else{//上面是0
if(sa==cur[y-1]){//上面连通块一样
//填0
add(dp[stat[msk][bl][wh]][sa],res);
//填1
int nb,nw;
color[y]=1,cur[y]=4;
restruct(nb,nw);
add(dp[stat[msk^(1<<y)][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}else{//上面连通块不一样
//填01
if(cnt!=1){
int nb,nw;
color[y]=1,cur[y]=4;
restruct(nb,nw);
add(dp[stat[msk^(1<<y)][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}
//填1
int nb,nw,dc=cur[y-1],df=cur[y];
rd(j,8)if(!color[j]&&cur[j]==df)cur[j]=dc;
restruct(nb,nw);
add(dp[stat[msk][nb][nw]][(y==8?cur[0]:cur[y+1])],res);
}
}
}
}
}
signed main(){
ios::sync_with_stdio(false);
cin.tie(0);
cin>>n>>m>>P;
init_div();
init_sta();
rd(i,(1<<m))f[stat[i][bel[cb[0][i]]-1][bel[cb[1][i]]-1]][0]=1;
rep(x,1,n-1)rep(y,0,m-2){
solve(y);
rd(i,stac)rd(j,4)f[i][j]=dp[i][j],dp[i][j]=0;
}
int ans=0;
rd(i,(1<<m))rd(j,4)add(ans,f[stat[i][0][0]][j]);
cout<<ans<<endl;
return 0;
}
//Crayan_r
便于实现的思路
既然一开始，我们就把状态设计的过于精简了，我们不如用别的方式。比如，我们可以用 \(Bell_8\)，\(Bell_8\) 其实也不大，只有 \(4000\) 左右，但它胜在从集合划分可以唯一映射到黑白染色方案。轮廓线 \(dp\) 的话，虽然会导致轮廓线处出现相同颜色不同集合，但是我们可以直接用一个 bool 变量记录在状态里。
然后，\(8\) 的集合划分也同样无序啊？但是我们可以把它编码了。首先，我们默认的集合划分是一个序列，表示每一位属于哪个集合。因为集合之间是不区分的，所以我们可以钦定第一个出现的是 \(0\)，第二个是 \(1\)。。。这样。我们发现，集合划分映射到一个长度 \(8\) 的串，第 \(i\) 位的值域是 \([0,i)\)。
这样，我们就可以哈希编码了。编码的总值域是 \(8!=40320\)，是上一个做法的 \(10\) 倍还要多。但是多出的这部分无关紧要，因为在常数部分就能补回来，而且稍微调大的状态数带来的是好写的转移。
然后，转移就比较方便。我们可以把转移拆成解码、转移和编码，解码部分就是把当前的状态解出黑白染色方案和集合划分方案，在同一列中，我们钦定第一个为黑色。转移的时候，通过枚举当前选择的颜色 \(t\)，存在四种情况：
-
和上面相同，和左面相同，则合并上面和左面的集合
-
和上面相同，和左面不同，首先要左面的集合在轮廓上还有其他出现，然后覆盖掉当前位置
-
和上面不同，和左面相同，等于什么都没有改变
-
和上面不同，和左面不同，还是要左边的集合在轮廓上有其他出现，然后建立新的划分集合。
这些转移都非常方便，可以两行左右写完。
然后着重是统计答案的问题。为了不算重，我们钦定凡是一列全是同一种颜色的都不再往下转移或者统计答案。也就是我们只统计整个图形黑白都有的列。
我们在每次 \(dp\) 完一列之后，先统计所有的黑白都有的答案，看看能否右边剩下的全填白色或者黑色。然后加入新的状态——左边全都是黑色，进行下一轮转移。尤其注意的是，在最后一列统计答案的时候，因为右边不会有东西了，所以我们强制要求有且仅有两个集合。其他的列上统计答案，只要某个颜色只有一个集合就可以计算一次贡献。
然后因为我们一开始钦定了左下角的颜色是黑色，所以方案乘以二。还要加上全黑或全白，得到最终答案。
ll n,m,P;
inline void add(ll &x,ll y){
x+=y;
if(x>=P)x-=P;
}
int a[9],x[8],cnt[16];
inline int incode(){
for(int i=0;i<16;i++)cnt[i]=0;
int tp=0;
rd(i,m){
if(!cnt[x[i]])cnt[x[i]]=++tp;
x[i]=cnt[x[i]]-1;
}
int res=0;
rep(i,0,m-1)res=res*(i+1)+x[i];
return res;
}
inline void decode(int c,int y){
per(i,0,m-1){
x[i]=c%(i+1);
c/=(i+1);
}
a[0]=0;
rep(i,1,m-1)if(x[i]!=x[i-1]){
a[i]=a[i-1]^(i!=y);
}else a[i]=a[i-1];
}
ll dp[50005][2],f[50005][2],N=50000;
inline void solve(int i){
rd(msk,N)rd(j,2)if(f[msk][j])rd(t,2){
decode(msk,i*j);
int cnt=0;
rd(p,m)if(x[p]==x[i])cnt++;
if(i!=0&&a[i-1]==t&&a[i]==t){
int fg=x[i-1],ft=x[i];
rd(p,m)if(x[p]==ft)x[p]=fg;
add(dp[incode()][a[i+1]==t],f[msk][j]);
}else if(i!=0&&a[i-1]==t&&a[i]!=t){
if(cnt!=1){
x[i]=x[i-1];
add(dp[incode()][a[i+1]==t],f[msk][j]);
}
}else if(a[i]==t){
add(dp[msk][a[i+1]==t],f[msk][j]);
}else if(cnt!=1){
x[i]=8;
add(dp[incode()][a[i+1]==t],f[msk][j]);
}
}
}
inline ll calc(){
ll res=0;
for(int cur=1;cur<N;cur++){
decode(cur,0);
int la=-1,lb=-1,fa=1,fb=1;
for(int i=0;i<m;i++){
if(a[i]){
if(la==-1)la=x[i];
else if(la!=x[i])fa=0;
}else{
if(lb==-1)lb=x[i];
else if(lb!=x[i])fb=0;
}
}
if(fa)add(res,f[cur][0]);
if(fb)add(res,f[cur][0]);
}
// cout<<res<<endl;
return res;
}
signed main(){
ios::sync_with_stdio(false);
cin.tie(0);
cin>>n>>m>>P;
if(n==1){
cout<<2*m%P<<endl;
return 0;
}
if(m==1){
cout<<2*n%P<<endl;
return 0;
}
for(int msk=0;msk<(1<<(m-1));msk++){
x[0]=0;
rd(i,m)a[i]=(msk>>i&1);
rep(i,1,m-1)x[i]=x[i-1]+(a[i]!=a[i-1]);
f[incode()][0]=1;
}a[m]=114514;
ll ans=calc();
rep(_,1,n-1){
rep(i,0,m-1){
solve(i);
x[0]=0,x[1]=1;
rep(i,0,N-1)rd(j,2)f[i][j]=dp[i][j],dp[i][j]=0;
}
if(_!=n-1)add(ans,calc());
f[0][0]=1;
}
for(int cur=1;cur<N;cur++){
decode(cur,0);
int la=-1,lb=-1,fa=1,fb=1;
for(int i=0;i<m;i++){
if(a[i]){
if(la==-1)la=x[i];
else if(la!=x[i])fa=0;
}else{
if(lb==-1)lb=x[i];
else if(lb!=x[i])fb=0;
}
}
if(fa&&fb&&f[cur][0]){
add(ans,f[cur][0]);
}
}
ans=(ans+n)%P;
cout<<(ans*2)%P<<endl;
return 0;
}
//Crayan_r