【NOIP2016提高A组五校联考4】square

Description

Input

Output

对于每个询问，输出答案

Sample Input

3 4
1 1 0 1
0 1 1 0
0 1 1 0
5
1 1 2 3
2 1 3 2
3 2 3 4
1 1 3 4
1 2 3 4

Sample Output

1
1
1
2
2

Data Constraint

Solution

多么良心的出题人！给出每个数据的大小，意思就是让你多水分。而且还注明了要卡常，而且还自己用C和P都打过了，确保PASCAL选手不因为语言问题被坑！多么好的出题人！
接着就恶心了。
一个显然的DP
设$f[i,j]$表示以(i,j)为右下角的正方形，边长最大多大
当(i,j)没树时，$f[i,j]=min(f[i-1,j-1],f[i,j-1],f[i-1,j])+1$
想到这就有很多分了
然后用一个$G[k,l,i,j]$表示以(i,j)为右下角，向左延伸$2^k$向上延伸$2^l$所表示的矩形的最大值，二维RMQ，就解决了
为什么要k,l在前呢？因为这样可以使访问内存时访问的位置尽可能集中从而减小常数。
加上读入优化之类的就可以切了
但是，如果你在求答案时发现了BUG：f[i,j]表示的正方形可能大过了给定的范围，怎么办？二分答案！

Code

#include<cstdio>
#include<algorithm>
#include<cmath>
#include<cstring>
#define fo(i,a,b) for(int i=a;i<=b;i++)
#define N 1020
using namespace std;
int n,m,f[N][N],g[10][10][N][N],a[N][N],e[15],lg[N];
int pd(int x1,int y1,int x2,int y2)
{
    int x=lg[x2-x1+1],y=lg[y2-y1+1];
    return max(max(g[x][y][x2][y2],g[x][y][x1+e[x]-1][y2]),max(g[x][y][x2][y1+e[y]-1],g[x][y][x1+e[x]-1][y1+e[y]-1]));
}
void read(int &n)
{
    char ch=' ';int q=0,w=1;
    for(;(ch!='-')&&((ch<'0')||(ch>'9'));ch=getchar());
    if(ch=='-') w=-1,ch=getchar();
    for(;ch>='0' && ch<='9';ch=getchar()) q=q*10+ch-48;n=q*w;
}
int main()
{
    freopen("square.in","r",stdin);freopen("square.out","w",stdout);
    e[0]=1;fo(i,1,14) e[i]=e[i-1]*2;
    fo(i,1,1001) lg[i]=(int)log2(i);
    read(n);read(m);
    fo(i,1,n)
    {
        fo(j,1,m)
        {
            read(a[i][j]);
            if(a[i][j]==1) f[i][j]=min(min(f[i-1][j-1],f[i-1][j]),f[i][j-1])+1,g[0][0][i][j]=f[i][j];
        }
    }
    fo(l,1,lg[m]) fo(i,1,n) fo(j,1,m) if(j-e[l-1]>=0) g[0][l][i][j]=max(g[0][l-1][i][j],g[0][l-1][i][j-e[l-1]]);
    fo(k,1,lg[n]) fo(i,1,n) if(i-e[k-1]>=0) fo(j,1,m) g[k][0][i][j]=max(g[k-1][0][i][j],g[k-1][0][i-e[k-1]][j]);    
    fo(k,1,lg[n]) fo(l,1,lg[m]) fo(i,1,n) if(i-e[k-1]>=0) fo(j,1,m) if(j-e[l-1]>=0) g[k][l][i][j]=max(max(g[k-1][l-1][i][j],g[k-1][l-1][i-e[k-1]][j-e[l-1]]),max(g[k-1][l-1][i-e[k-1]][j],g[k-1][l-1][i][j-e[l-1]]));
    int ac;read(ac);
    for(;ac;ac--)
    {
        int x1,x2,y1,y2;read(x1);read(y1);read(x2);read(y2);
        int l=0,r=min(x2-x1+1,y2-y1+1);
        for(;l+1<r;)
        {
            int md=(l+r)/2;int x3=md+x1-1,y3=md+y1-1;
            if(pd(x3,y3,x2,y2)>=md) l=md;else r=md;
        }
        int x3=r+x1-1,y3=r+y1-1;
        if(pd(x3,y3,x2,y2)>=r) l=r;
        printf("%d\n",l);
    }
}