https://www.51nod.com/onlineJudge/questionCode.html#!problemId=1303
所求=总矩阵对数-相离矩阵对数-相包含矩阵对数
总矩阵对数=${矩阵数 \choose 2}$
考虑如何求有多少个全1矩阵。
先枚举右下角的端点$(i,j)$，维护以$(i,j)$为右下角的全1矩阵个数
枚举行$i$，从左到右枚举列$j$
用单调栈维护$k=1..j$到$j$的$up$的$min$ ($up[i][j]$表示点$(i,j)$向上有连续多少个1)
以$(i,j)$为右下角的矩阵个数=$min_k$的和，可以维护。
考虑如何求相包含矩阵对数。
枚举较大的矩阵($n*m$)，则较小矩阵有${n+1 \choose 2} * {m+1 \choose 2}$个
所以用类似上面的方法，维护$1*f(m),n*f(m),n^2*f(m)$的和即可。($f(m)为{m+1 \choose 2} 的前缀和$)
考虑如何求相离矩阵对数。
即左右相离+上下相离-左上右下相离-右上左下相离
求出每个点作为四个方向的端点的矩阵个数，之后就轻松解决了。
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define pb push_back
#define rep(i,l,r) for(int i=l;i<=r;++i)
#define per(i,r,l) for(int i=r;i>=l;--i)
const int N=1e3+5,D=1e9+7;
ll mi(ll x,int y=D-2)
{
ll ans=1;
while(y)
{
if(y&1)ans=ans*x%D;
x=x*x%D;y>>=1;
}
return ans;
}
int n,m;
char b[N][N];
int up[N][N],down[N][N];
ll S1[N],S2[N],f[N],f2[N];
//S1[n]=sigma i; S2[n]=sigma i^2; f[n]=sigma S[i]
ll a[N],h[N];
ll cnt[N][N],cnt1[N][N],cnt2[N][N],cnt3[N][N];
//cnt 右下角 ;cnt1 左下角 ;cnt2 右上角 ;cnt3 左上角
ll c[N],c1[N];
int main()
{
freopen("1.in","r",stdin);//freopen("2.out","w",stdout);
cin>>n>>m;
rep(i,1,n)scanf("%s",b[i]+1);
int mn=max(n,m);
rep(i,1,mn){S1[i]=S1[i-1]+i;f[i]=f[i-1]+S1[i];}
rep(i,1,mn){S2[i]=S2[i-1]+i*i;(f2[i]=f2[i-1]+S2[i])%D;}
rep(i,1,n)
rep(j,1,m)
if(b[i][j]=='1')up[i][j]=up[i-1][j]+1;
per(i,n,1)
rep(j,1,m)
if(b[i][j]=='1')down[i][j]=down[i+1][j]+1;
ll ans=0;
//+C2(tot_cnt)
//-相包含
//-横相离-纵相离+左上右下相离+右上左下相离
//-相包含
rep(i,1,n)
{
int top=0,ncnt=0;
//ncnt sigma h_k
ll s0=0,s1=0,s2=0;
//s0=sigma f[h_k]
//s1=sigma (j-k+1)*f[h_k]
//s2=sigma (j-k+1)^2*f[h_k]
rep(j,1,m)
{
(s2+=2*s1+s0)%=D;
(s1+=s0)%=D;
while(h[top]>up[i][j])
{
ncnt-=h[top]*(a[top]-a[top-1]);
(s0-=f[h[top]]*(a[top]-a[top-1]))%=D;
(s1-=f[h[top]]*(S1[j-a[top-1]]-S1[j-a[top]]))%=D;
(s2-=f[h[top]]*(S2[j-a[top-1]]-S2[j-a[top]]))%=D;
--top;
}
h[++top]=up[i][j];
a[top]=j;
ncnt+=h[top]*(a[top]-a[top-1]);
(s0+=f[h[top]]*(a[top]-a[top-1]))%=D;
(s1+=f[h[top]]*(S1[j-a[top-1]]-S1[j-a[top]]))%=D;
(s2+=f[h[top]]*(S2[j-a[top-1]]-S2[j-a[top]]))%=D;
cnt[i][j]=ncnt;
(ans+=s2+s1)%=D;
}
}
(ans*=mi(2))%=D;
ans=-ans;
//get cnt2
rep(i,1,n)
{
int top=0,ncnt=0;
rep(j,1,m)
{
while(h[top]>down[i][j])
{
ncnt-=h[top]*(a[top]-a[top-1]);
--top;
}
h[++top]=down[i][j];
a[top]=j;
ncnt+=h[top]*(a[top]-a[top-1]);
cnt2[i][j]=ncnt;
}
}
a[0]=m+1;
//get cnt1
rep(i,1,n)
{
int top=0,ncnt=0;
per(j,m,1)
{
while(h[top]>up[i][j])
{
ncnt-=h[top]*(a[top-1]-a[top]);
--top;
}
h[++top]=up[i][j];
a[top]=j;
ncnt+=h[top]*(a[top-1]-a[top]);
cnt1[i][j]=ncnt;
}
}
//get cnt3
rep(i,1,n)
{
int top=0,ncnt=0;
per(j,m,1)
{
while(h[top]>down[i][j])
{
ncnt-=h[top]*(a[top-1]-a[top]);
--top;
}
h[++top]=down[i][j];
a[top]=j;
ncnt+=h[top]*(a[top-1]-a[top]);
cnt3[i][j]=ncnt;
}
}
//-横相离
rep(i,1,n)
rep(j,1,m) c[j]+=cnt[i][j];
rep(i,1,n)
rep(j,1,m) c1[j]+=cnt1[i][j];
ll sum=0;
rep(i,1,m)sum+=c[i];
sum%=D;
//----- +C2(tot_cnt)
(ans+=(sum+1)*sum/2)%=D;
//-----
per(i,m,1)
{
(sum-=c[i])%=D;
(ans-=sum*c1[i])%=D;
}
//-纵相离
memset(c,0,sizeof(c));memset(c1,0,sizeof(c1));
rep(i,1,n)
rep(j,1,m) c[i]+=cnt[i][j];
rep(i,1,n)
rep(j,1,m) c1[i]+=cnt2[i][j];
sum=0;
rep(i,1,n)sum+=c[i];
sum%=D;
per(i,n,1)
{
(sum-=c[i])%=D;
(ans-=sum*c1[i])%=D;
}
//+左上右下相离
rep(i,1,n)
rep(j,1,m)cnt[i][j]+=cnt[i][j-1];
rep(i,1,n)
rep(j,1,m)(cnt[i][j]+=cnt[i-1][j])%=D;
rep(i,1,n)
rep(j,1,m)(ans+=cnt3[i][j]*cnt[i-1][j-1])%=D;
//+左下右上相离
rep(i,1,n)
per(j,m,1)cnt1[i][j]+=cnt1[i][j+1];
rep(i,1,n)
rep(j,1,m)(cnt1[i][j]+=cnt1[i-1][j])%=D;
rep(i,1,n)
rep(j,1,m)(ans+=cnt2[i][j]*cnt1[i-1][j+1])%=D;
cout<<(ans%D+D)%D;
}
鄂公网安备 42010202000505 号