数学（容斥计数）：LNOI 2016 方

Description

上帝说，不要圆，要方，于是便有了这道题。由于我们应该方，而且最好能够尽量方，所以上帝派我们来找正方形

上帝把我们派到了一个有N行M列的方格图上，图上一共有(N+1)×(M+1)个格点，我们需要做的就是找出这些格点形

成了多少个正方形（换句话说，正方形的四个顶点都是格点）。但是这个问题对于我们来说太难了，因为点数太多

了，所以上帝删掉了这(N+1)×(M+1)中的K个点。既然点变少了，问题也就变简单了，那么这个时候这些格点组成

了多少个正方形呢？

Input

第一行三个整数 N, M, K， 代表棋盘的行数、 列数和不能选取的顶点个数。 保证 N, M >= 1， K <=(N + 1) ×

(M + 1)。约定每行的格点从上到下依次用整数 0 到 N 编号，每列的格点依次用 0到 M 编号。接下来 K 行，每

行两个整数 x,y 代表第 x 行第 y 列的格点被删掉了。保证 0 <=x <=N<=10^6, 0 <=y<=M<=10^6，K<=2*1000且不

会出现重复的格点。

Output

 仅一行一个正整数， 代表正方形个数对 100000007（ 10^8 + 7） 取模之后的值

Sample Input

2 2 4
1 0
1 2
0 1
2 1

Sample Output

1

　　这道题因为k很小可以用容斥， 所有正方形数-Σ每个点（在正方形上）的正方形数+Σ每个无序点对（在正方形上）的正方形数-Σ每个无序三个元组（在正方形上）的正方形数+Σ每个四元组（在正方形上）的正方形数

　　然后我打的复杂化了，不过很好懂，几乎是望文生义……

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdio>
using namespace std;
const int N=2010;
const int Mod=100000007;
typedef long long LL;
LL Abs(LL x){return x<0?-x:x;}
LL n,m,k,tot,x[N],y[N];
LL Ask(LL a,LL b,LL c){
    LL ret,L;
    if(!a||!b||!c)return 0;
    if(a>c-1)a=c-1;
    if(b>c-1)b=c-1;
    ret=(2*c-a-1)*a/2;
    if(c-1>=b+1){
        L=min(a,c-b-1);
        ret+=b*L;
        ret-=(2*c-L-1)*L/2;
    }
    return (ret%Mod+Mod)%Mod;
}
LL Solve1(){
    LL ret=0,L,R,U,D;
    for(int i=1;i<=k;i++){
        L=min(x[i]-1,y[i]-1);ret+=L;
        L=min(n+1-x[i],y[i]-1);ret+=L;
        L=min(x[i]-1,m+1-y[i]);ret+=L;
        L=min(n+1-x[i],m+1-y[i]);ret+=L;
    }
    for(int i=1;i<=k;i++){
        U=x[i]-1;D=n+1-x[i];
        L=y[i]-1;R=m+1-y[i];
        ret+=Ask(L,R,U);
        ret+=Ask(L,R,D);
        ret+=Ask(U,D,L);
        ret+=Ask(U,D,R);
    }
    return ret%Mod;
}


LL Solve2(){
    LL ret=0,L;
    LL px,py,dx,dy;
    for(int i=1;i<=k;i++)
        for(int j=i+1;j<=k;j++){
            dx=x[j]-x[i];
            dy=y[j]-y[i];
            
            px=x[j]+dy;py=y[j]-dx;
            if(px>=1&&px<=n+1&&py>=1&&py<=m+1){
                px=px-dx;py=py-dy;
                if(px>=1&&px<=n+1&&py>=1&&py<=m+1)ret+=1;    
            }
            
            px=x[j]-dy;py=y[j]+dx;
            if(px>=1&&px<=n+1&&py>=1&&py<=m+1){
                px=px-dx;py=py-dy;
                if(px>=1&&px<=n+1&&py>=1&&py<=m+1)ret+=1;    
            }
            
            if((x[i]+x[j])%2==(y[i]+y[j])%2){
                dx=2*x[i]-(x[i]+x[j]);
                dy=2*y[i]-(y[i]+y[j]);
                px=(x[i]+x[j]+dy)/2;
                py=(y[i]+y[j]-dx)/2;
                if(px>=1&&px<=n+1&&py>=1&&py<=m+1){
                    dx=-dx;dy=-dy;
                    px=(x[i]+x[j]+dy)/2;
                    py=(y[i]+y[j]-dx)/2;
                    if(px>=1&&px<=n+1&&py>=1&&py<=m+1)ret+=1;
                }
            }
        }
    return ret%Mod;        
}

LL hshx[N],hshy[N];
int cntx,cnty,vis[N][N];

void Prepare(){
    for(int i=1;i<=k;i++){
        hshx[++cntx]=x[i];
        hshy[++cnty]=y[i];
    }
    sort(hshx+1,hshx+cntx+1);
    cntx=unique(hshx+1,hshx+cntx+1)-hshx-1;
    sort(hshy+1,hshy+cnty+1);
    cnty=unique(hshy+1,hshy+cnty+1)-hshy-1;
    for(int i=1;i<=k;i++){
        LL x_=lower_bound(hshx+1,hshx+cntx+1,x[i])-hshx;
        LL y_=lower_bound(hshy+1,hshy+cnty+1,y[i])-hshy;
        vis[x_][y_]=i;
    }
}

int Q(LL x,LL y){
    LL x_=lower_bound(hshx+1,hshx+cntx+1,x)-hshx;
    LL y_=lower_bound(hshy+1,hshy+cnty+1,y)-hshy;
    return (hshx[x_]==x&&hshy[y_]==y)*vis[x_][y_];
}

LL Solve3(){
    LL ret=0;
    LL px1,py1,px2,py2,dx,dy;
    for(int i=1;i<=k;i++)
        for(int j=i+1;j<=k;j++){
            dx=x[j]-x[i];
            dy=y[j]-y[i];
            
            px1=x[j]+dy;py1=y[j]-dx;
            if(px1>=1&&px1<=n+1&&py1>=1&&py1<=m+1){
                px2=px1-dx;py2=py1-dy;
                if(px2>=1&&px2<=n+1&&py2>=1&&py2<=m+1){
                    if(Q(px1,py1)>j)ret+=1;
                    if(Q(px2,py2)>j)ret+=1;
                }
            }
            
            px1=x[j]-dy;py1=y[j]+dx;
            if(px1>=1&&px1<=n+1&&py1>=1&&py1<=m+1){
                px2=px1-dx;py2=py1-dy;
                if(px2>=1&&px2<=n+1&&py2>=1&&py2<=m+1){
                    if(Q(px1,py1)>j)ret+=1;
                    if(Q(px2,py2)>j)ret+=1;
                }
            }
            
            if((x[i]+x[j])%2==(y[i]+y[j])%2){
                dx=2*x[i]-(x[i]+x[j]);
                dy=2*y[i]-(y[i]+y[j]);
                px1=(x[i]+x[j]+dy)/2;
                py1=(y[i]+y[j]-dx)/2;
                if(px1>=1&&px1<=n+1&&py1>=1&&py1<=m+1){
                    dx=-dx;dy=-dy;
                    px2=(x[i]+x[j]+dy)/2;
                    py2=(y[i]+y[j]-dx)/2;
                    if(px2>=1&&px2<=n+1&&py2>=1&&py2<=m+1){
                        if(Q(px1,py1)>j)ret+=1;
                        if(Q(px2,py2)>j)ret+=1;
                    }
                }
            }
        }
    return ret%Mod;
}

LL Solve4(){
    LL ret=0;
    LL px1,py1,px2,py2,dx,dy;
    for(int i=1;i<=k;i++)
        for(int j=i+1;j<=k;j++){
            dx=x[j]-x[i];
            dy=y[j]-y[i];
            
            px1=x[j]+dy;py1=y[j]-dx;
            if(px1>=1&&px1<=n+1&&py1>=1&&py1<=m+1){
                px2=px1-dx;py2=py1-dy;
                if(px2>=1&&px2<=n+1&&py2>=1&&py2<=m+1)
                    if(Q(px1,py1)>j&&Q(px2,py2)>j)ret+=1;
            }
            
            px1=x[j]-dy;py1=y[j]+dx;
            if(px1>=1&&px1<=n+1&&py1>=1&&py1<=m+1){
                px2=px1-dx;py2=py1-dy;
                if(px2>=1&&px2<=n+1&&py2>=1&&py2<=m+1)
                    if(Q(px1,py1)>j&&Q(px2,py2)>j)ret+=1;
            }
            
            if((x[i]+x[j])%2==(y[i]+y[j])%2){
                dx=2*x[i]-(x[i]+x[j]);
                dy=2*y[i]-(y[i]+y[j]);
                px1=(x[i]+x[j]+dy)/2;
                py1=(y[i]+y[j]-dx)/2;
                if(px1>=1&&px1<=n+1&&py1>=1&&py1<=m+1){
                    dx=-dx;dy=-dy;
                    px2=(x[i]+x[j]+dy)/2;
                    py2=(y[i]+y[j]-dx)/2;
                    if(px2>=1&&px2<=n+1&&py2>=1&&py2<=m+1)
                        if(Q(px1,py1)>j&&Q(px2,py2)>j)ret+=1;
                }
            }
        }
    return ret%Mod;
}

int main(){
    freopen("square.in","r",stdin);
    freopen("square.out","w",stdout);
    scanf("%lld%lld%lld",&n,&m,&k);
    for(LL L=min(m,n);L>=1;L--){
        tot+=(n-L+1)*(m-L+1)%Mod*L%Mod;
        if(tot>=Mod)tot-=Mod;
    }
    for(int i=1;i<=k;i++){
        scanf("%lld%lld",&x[i],&y[i]);
        x[i]+=1;y[i]+=1;
    }
    Prepare();
    tot-=Solve1();tot%=Mod;
    tot+=Solve2();tot%=Mod;
    tot-=Solve3();tot%=Mod;
    tot+=Solve4();tot%=Mod;    
    tot=(tot+Mod)%Mod;
    printf("%lld
",tot);
    return 0;
}

　　运算有些多，常数比较大，不开O2只有50分，优化的话，可以把LL改成int，Mod的时候有些地方可以改成减。