BZOJ1227 SDOI2009 虔诚的墓主人 树状数组+组合数学

题意：给定一个N*M的01矩阵，设l_i,j,r_i,j,u_i,j,d_i,j分别为(i,j)正上，正下，正左，正右1的数量，求$\sum\limits_{(i,j) \equiv 0} {C_{{l_{i,j}}}^KC_{{r_{i,j}}}^KC_{{u_{i,j}}}^KC_{{d_{i,j}}}^K} $，其中1的数量≤100000

题解：首先离散化，设N M为离散化后的行列数。先把l r u d全算出来，由于对于固定的一行i，有\[\sum\limits_{j = 1}^M {C_{{l_{i,j}}}^KC_{{r_{i,j}}}^KC_{{u_{i,j}}}^KC_{{d_{i,j}}}^K}  = \sum\limits_{(i,a) = 1,(i,b) = 1,b > a} {\left( {C_{{l_{i,a}}}^KC_{{r_{i,b}}}^K\sum\limits_{j = a + 1}^{b - 1} {C_{{u_{i,j}}}^KC_{{d_{i,j}}}^K} } \right)} \]

因此我们可以枚举每一个1，问题的关键集中在如何求第二个求和。设(x,i)=1,(y,i)=1且y为x正上方的第一个1，那么有\[\begin{array}{l}
{u_{y,i}} = {u_{x,i}} - 1\\
{d_{y,i}} = {d_{x,i}} + 1
\end{array}\]

我们用树状数组来维护求和，没枚举到一个1，就删除其上第一个点的求和，加入当前点的求和。由于取余的数很特殊，所以我们可以自然溢出，最后&((2^31)-1)即可。

尼玛组合数写炸了WA了一发……死活调不出来

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <iostream>
#include <algorithm>
using namespace std;
#define P 2147483647
#define lowbit(x) (x&(-x))

const int MAXK=10+2;
const int MAXW=100000+2;
struct Node{ int x,y,m;}t[MAXW];
int N,M,W,K;
int s[2][MAXW],l[MAXW],r[MAXW],u[MAXW],d[MAXW];
int a[MAXW],Ans,Com[MAXW][MAXK];

bool cmp1(Node a,Node b){
    if(a.x==b.x) return a.y<b.y;
    return a.x<b.x;
}

bool cmp2(Node a,Node b){
    if(a.y==b.y) return a.x<b.x;
    return a.y<b.y;
}

void Update(int p,int x){
    while(p<=M) a[p]+=x,p+=lowbit(p);
}

int Query(int p){
    int Ret=0;
    while(p) Ret+=a[p],p-=lowbit(p);
    return Ret;
}

int C(int n,int m){ return Com[n][m];}

int main(){
    scanf("%d %d %d",&N,&M,&W);
    memset(t,0,sizeof(t));
    for(int i=1;i<=W;i++) scanf("%d %d",&t[i].x,&t[i].y);
    scanf("%d",&K);

    sort(t+1,t+W+1,cmp1),N=0;
    for(int i=1;i<=W;i++) t[i].m=i;
    for(int i=1,x=-1;i<=W;i++){
        if(x!=t[i].x) ++N;
        x=t[i].x,t[i].x=N,s[0][N]++;
    }
    for(int i=1,c;i<=W;i++){
        if(t[i].x!=t[i-1].x) c=-1;
        c++,l[i]=c,r[i]=s[0][t[i].x]-c-1;
    }
    sort(t+1,t+W+1,cmp2),M=0;
    for(int i=1,y=-1;i<=W;i++){
        if(y!=t[i].y) ++M;
        y=t[i].y,t[i].y=M,s[1][M]++;
    }
    for(int i=1,c;i<=W;i++){
        if(t[i].y!=t[i-1].y) c=-1;
        c++,u[t[i].m]=c,d[t[i].m]=s[1][t[i].y]-c-1;
    }

    for(int i=0,j;i<=W;i++)
        for(j=1,Com[i][0]=1;j<=min(i,K);j++)
            Com[i][j]=Com[i-1][j-1]+Com[i-1][j];

    sort(t+1,t+W+1,cmp1);
    for(int i=1;i<=W;i++){
        Update(t[i].y,C(u[i]+1,K)*C(d[i],K)-C(u[i],K)*C(d[i]+1,K));
        while(i<W && t[i].x==t[i+1].x){
            i++;
            Ans+=C(l[i],K)*C(r[i-1],K)*(Query(t[i].y-1)-Query(t[i-1].y));
            Update(t[i].y,C(u[i]+1,K)*C(d[i],K)-C(u[i],K)*C(d[i]+1,K));
        }
    }
    printf("%d\n",Ans&P);

    return 0;
}