题解：P14476 矩阵绘画

· · 题解

如果对合法染色方式计数，那么是容易的。考虑不相交的限制，提取出激活的 a 的前缀 \max，对这些前缀 \max 做 dp，此时就将矩阵划分成了一些区域，每个区域要么全选列要么全选行，这样是容易计数的，处理 [1,i] 中所有后缀中 b_j<a_i 的信息即可，再预处理一下每个行，列的 0 和 1 信息，可以 \mathcal O(n^2) 的复杂度算。

问题在于重复计数。考虑一个重复算的区域，要么是前缀全选行，要么是前缀全选列（一个 2\times 2 内部全为 1 的区域，考虑两行均被行染色以及两行均被列染色）。a_1=0 时压根不会数重，进而，猜测每个矩阵只会被至多两种方案算到。这样我们容斥：用总操作方案数，减去会被算到两次的矩阵个数，就是最终答案。

进而，只关心左上角 x\times y 的矩阵内部，要求 \max_{i\leq x}a_i\leq y,\max_{i\leq y}b_i\leq x，并且全选行与全选列构成的 01 矩阵完全相同。然后再乘上剩余部分的染色方案数即可，一个必要条件是 a,b 单调不增，进而发现，这个区域只能取到 x=b_1,y=a_1。这里再做 dp 可以沿用对 a 的前缀 \max 计数的思路，但是要修改 a_1 之后第一个前缀 \max 的值。

时间复杂度 \mathcal O(n^2)。

#include<bits/stdc++.h>
#define ull unsigned long long
#define ll long long
#define uint unsigned
#define pb push_back
#define mkp make_pair
#define fi first
#define se second
#define inf 1000000000
#define infll 1000000000000000000ll
#define pii pair<int,int>
#define rep(i,a,b,c) for(int i=(a);i<=(b);i+=(c))
#define per(i,a,b,c) for(int i=(a);i>=(b);i-=(c))
#define F(i,a,b) for(int i=a,i##end=b;i<=i##end;i++)
#define dF(i,a,b) for(int i=a,i##end=b;i>=i##end;i--)
#define eb emplace_back
#define SZ(x) ((int)x.size())
#define all(x) x.begin(),x.end()
using namespace std;
bool ST;
void fre(){
    // freopen(".in","r",stdin),freopen(".out","w",stdout);
}
inline int lowbit(int x){ return x&(-x); }
template<typename T>inline void chkmax(T &x,const T &y){ x=std::max(x,y); }
template<typename T>inline void chkmin(T &x,const T &y){ x=std::min(x,y); }
const int mod=1e9+7,maxn=5005;
inline int qpow(int x,ll y){ int res=1; for(;y;y>>=1,x=1ll*x*x%mod)if(y&1)res=1ll*res*x%mod; return res; }
inline void inc(int &x,const int y){ x=(x+y>=mod)?(x+y-mod):(x+y); }
inline void dec(int &x,const int y){ x=(x>=y)?(x-y):(x+mod-y); }
inline int add(const int x,const int y){ return (x+y>=mod)?(x+y-mod):(x+y); }
inline int sub(const int x,const int y){ return (x>=y)?(x-y):(x+mod-y); }
string str[maxn];
int n,a[maxn],b[maxn],p0[maxn][maxn],p1[maxn][maxn],h[maxn][maxn],dp[maxn],g[maxn],f[maxn];
inline int q0(int xl,int xr,int yl,int yr){ return (xl>xr||yl>yr)?0:p0[xr][yr]-p0[xl-1][yr]-p0[xr][yl-1]+p0[xl-1][yl-1]; }
inline int q1(int xl,int xr,int yl,int yr){ return (xl>xr||yl>yr)?0:p1[xr][yr]-p1[xl-1][yr]-p1[xr][yl-1]+p1[xl-1][yl-1]; }
int sol1(){
    f[0]=1,a[n+1]=n;
    F(i,1,n)if(a[i]){
        h[i][i]=1;
        F(j,i+1,n+1){
            h[i][j]=0;
            if(a[j]>=a[i]){
                if(!q1(j,j,1,a[i]))h[i][j]=h[i][j-1];
                continue;
            }
            if(q1(j,j,a[j]+1,a[i])>0)continue;
            if(q0(j,j,1,a[j])==0)inc(h[i][j],h[i][j-1]);
            if(a[j]&&q1(j,j,1,a[j])==0)inc(h[i][j],h[i][j-1]);
        }
    }
    F(i,1,n+1)if(a[i]){
        g[a[i]+1]=1;
        dF(j,a[i],1){
            g[j]=0;
            if(b[j]>=i){
                if(!q1(1,i-1,j,j))g[j]=g[j+1];
                continue;
            }
            if(q1(b[j]+1,i-1,j,j)>0)continue;
            if(q0(1,b[j],j,j)==0)inc(g[j],g[j+1]);
            if(b[j]&&q1(1,b[j],j,j)==0)inc(g[j],g[j+1]);
        }
        F(j,0,i-1)if(a[j]<=a[i])inc(f[i],1ll*(j==0?1:h[j][i-1])*g[a[j]+1]%mod*f[j]%mod);
    }
    return f[n+1];
}
mt19937_64 eng(time(0));
uint rd[maxn][maxn];
bool chk(){
    if(!a[1]||!b[1])return 0;
    uint x=0,y=0;
    F(i,1,n)F(j,1,n)rd[i][j]=eng();
    F(i,1,a[1])F(j,1,b[i])x+=rd[j][i];
    F(i,1,b[1])F(j,1,a[i])y+=rd[i][j];
    return x==y;
}
int sol2(){
    memset(f,0,sizeof f),f[b[1]]=1;
    a[b[1]]=a[1];
    F(i,b[1],n)if(a[i]){
        h[i][i]=1;
        F(j,i+1,n+1){
            h[i][j]=0;
            if(a[j]>=a[i]){
                if(!q1(j,j,1,a[i]))h[i][j]=h[i][j-1];
                continue;
            }
            if(q1(j,j,a[j]+1,a[i])>0)continue;
            if(q0(j,j,1,a[j])==0)inc(h[i][j],h[i][j-1]);
            if(a[j]&&q1(j,j,1,a[j])==0)inc(h[i][j],h[i][j-1]);
        }
    }
    F(i,b[1]+1,n+1)if(a[i]){
        g[a[i]+1]=1;
        dF(j,a[i],1){
            g[j]=0;
            if(b[j]>=i){
                if(!q1(1,i-1,j,j))g[j]=g[j+1];
                continue;
            }
            if(q1(b[j]+1,i-1,j,j)>0)continue;
            if(q0(1,b[j],j,j)==0)inc(g[j],g[j+1]);
            if(b[j]&&q1(1,b[j],j,j)==0)inc(g[j],g[j+1]);
        }
        F(j,0,i-1)if(a[j]<=a[i])inc(f[i],1ll*(j==0?1:h[j][i-1])*g[a[j]+1]%mod*f[j]%mod);
    }
    return f[n+1];
}
void solve(){
    cin>>n;
    F(i,1,n)cin>>a[i];
    F(i,1,n)cin>>b[i];
    F(i,1,n)cin>>str[i],str[i]=' '+str[i];
    F(i,1,n)F(j,1,n)p0[i][j]=p0[i-1][j]+p0[i][j-1]-p0[i-1][j-1]+(str[i][j]=='0');
    F(i,1,n)F(j,1,n)p1[i][j]=p1[i-1][j]+p1[i][j-1]-p1[i-1][j-1]+(str[i][j]=='1');
    int ans=sol1();
    if(chk())dec(ans,sol2());
    cout<<ans<<'\n';
}
bool ED;
signed main(){
    fre(),ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
    int wzq=1;
    F(____,1,wzq)solve();
    cerr<<"time used: "<<(double)clock()/CLOCKS_PER_SEC<<endl;
    cerr<<"memory used: "<<abs(&ST-&ED)/1024.0/1024.0<<" MB"<<endl;
}
// g++ P14476.cpp -o a -std=c++14 -O2