bzoj1227 [SDOI2009]虔诚的墓主人（树状数组）

题目描述

小W是一片新造公墓的管理人。公墓可以看成一块N×M的矩形，矩形的每个格点，要么种着一棵常青树，要么是一块还没有归属的墓地。
当地的居民都是非常虔诚的基督徒，他们愿意提前为自己找一块合适墓地。为了体现自己对主的真诚，他们希望自己的墓地拥有着较高的虔诚度。
一块墓地的虔诚度是指以这块墓地为中心的十字架的数目。一个十字架可以看成中间是墓地，墓地的正上、正下、正左、正右都有恰好k棵常青树。
小W希望知道他所管理的这片公墓中所有墓地的虔诚度总和是多少。

输入文件

输入文件religious.in的第一行包含两个用空格分隔的正整数N和M，表示公墓的宽和长，因此这个矩形公墓共有(N+1) ×(M+1)个格点，左下角的坐标为(0, 0)，右上角的坐标为(N, M)。
第二行包含一个正整数W，表示公墓中常青树的个数。
第三行起共W行，每行包含两个用空格分隔的非负整数xi和yi，表示一棵常青树的坐标。输入保证没有两棵常青树拥有相同的坐标。
最后一行包含一个正整数k，意义如题目所示。

输出文件

仅包含一个非负整数，表示这片公墓中所有墓地的虔诚度总和。为了方便起见，答案对2,147,483,648取模。

样例输入

5 6
13
0 2
0 3
1 2
1 3
2 0
2 1
2 4
2 5
2 6
3 2
3 3
4 3
5 2
2

样例输出

6

样例说明

图中，以墓地(2, 2)和(2, 3)为中心的十字架各有3个，即它们的虔诚度均为3。其他墓地的虔诚度为0。

数据规模和约定

对于30%的数据，满足1 ≤ N, M ≤ 1,000。
对于60%的数据，满足1 ≤ N, M ≤ 1,000,000。
对于100%的数据，满足1 ≤ N, M ≤ 1,000,000,000，0 ≤ xi ≤ N，0 ≤ yi ≤ M，1 ≤ W ≤ 100,000，1 ≤ k ≤ 10。
存在50%的数据，满足1 ≤ k ≤ 2。
存在25%的数据，满足1 ≤ W ≤ 10000。

分析：
显然只有有点的一行（或一列）才有可能产生有贡献的墓地
因此只有W*W个坐标有贡献，所以我们先离散一下坐标

设一个点上下左右的常青树的个数分别是：u[i]，d[i]，l[i]，r[i]（可以直接预处理出来）
那么这个点的贡献就是：C(u[i] , K) * C(r[i] , K) * C(d[i] , K) * C(r[i] , K)

我们按照y坐标排序
这样就可以一行一行处理贡献了
在同一行中，两个相邻的常青树之间（a和b）的所有墓地都是有可能产生贡献
则ans+=C(l[a]+1(包括a)，k) * C(r[b]+1,k) * [乘上ab间的每一个点的C(u[i] , k) * C(d[i] , k)]
这样的复杂度是O(n^2)

Can we do better？
我们观察ab之间的C(u[i] , k) * C(d[i] , k)是一个区间和的形式
我们可以前缀和优化，
但是我们还要支持快速修改：树状数组
用树状数组维护a到b所有点的C(u[i],k)*C(d[i],k)之和

我们一行一行的处理，每次计算同一行中相邻结点之间的贡献
但是这里有一个小问题：
比如某列某行有一个点
在扫描这行之下的时候，这个点是算在u[i]里的，但是扫描这行之上时算在了d[i]中

于是我们从左往右处理某行的某一个点时，要将树状数组中该点横坐标位置上的数清除

#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#define ll long long

using namespace std;

const int p=2147483647;                 //2^31-1
const int N=100005;
int n,m,sum[N];
int K,W,X[N],Y[N],C[N][11],tree[N];
struct node{
    int x,y,u,d,l,r;
};
node a[N];

int cmp(const node &a,const node &b)   //按照y排序 
{
    if (a.y!=b.y) return a.y<b.y;
    else return a.x<b.x;
}

void pre()
{
    C[1][0]=1;C[1][1]=1;                         //预处理组合数 
    for (int i=2;i<=W;i++){
        C[i][0]=1; 
        for (int j=1;j<=K&&j<=i;j++) C[i][j]=C[i-1][j]+C[i-1][j-1];
    }

    sort(a+1,a+1+W,cmp);               //预处理每个点的u,d,l,r(包括本身) 
    int tmp=0;
    for (int i=1;i<=W;i++)
    {
        if (a[i].y==a[i-1].y) tmp++;
        else tmp=1;
        a[i].l=tmp; 
        sum[a[i].x]++; a[i].u=sum[a[i].x];
    }
    tmp=0;
    for (int i=W;i>=1;i--)
    {
        if (a[i].y==a[i+1].y) tmp++;
        else tmp=1;
        a[i].r=tmp;
        a[i].d=sum[a[i].x]-a[i].u;
    }
} 

void add(int x,int z)
{
    for (int i=x;i<=W;i+=(i&(-i))) tree[i]+=z;
}

int ask(int x)
{
    int ans=0;
    for (int i=x;i>0;i-=(i&(-i))) ans+=tree[i];
    return ans;
}

void solve()
{
    int ans=0;
    for (int i=1;i<=W;i++)
    {
        int q=ask(a[i].x)-ask(a[i].x-1);      //x位置原来的数 
        add(a[i].x,-q);
        add(a[i].x,C[a[i].u][K]*C[a[i].d][K]);
        if (a[i].y==a[i-1].y)                //一行 
            ans=ans+C[a[i-1].l][K]*C[a[i].r][K]*(ask(a[i].x-1)-ask(a[i-1].x));
    }
    printf("%d",ans&p);
}

int main()
{
    scanf("%lld%lld",&n,&m);
    scanf("%d",&W);
    for (int i=1;i<=W;i++)
        scanf("%d%d",&a[i].x,&a[i].y),X[i]=a[i].x,Y[i]=a[i].y;

    sort(X+1,X+1+W); sort(Y+1,Y+1+W);          //离散 
    for (int i=1;i<=W;i++)
        a[i].x=lower_bound(X+1,X+1+W,a[i].x)-X,
        a[i].y=lower_bound(Y+1,Y+1+W,a[i].y)-Y;

    scanf("%d",&K);
    pre();
    solve();

    return 0;
}