「SDOI2009」虔诚的墓主人 - 组合计数+树状数组+扫描

题目大意

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

---

题目分析

首先，对棋盘离散化。
通过递推预处理出$Left[],Right[],Up[],Down[]$分别表示每个常青树左右上下的常青树个数（不含本身）。
考虑矩阵的一行$y$中两个相邻的常青树$p,q(p\lt q)$，记$U[x,y]$表示墓地$(x,y)$上方的总常青树个数，$D[l,r]$表示下方（不含本身）。
那么这一段对答案的贡献即为：

$$C_{l[p]+1}^kC_{r[q]+1}^k\sum_{x=x_p+1}^{x_q-1}C_{U[x,y_p]}^kC_{D[x,y_p]}^k$$

考虑沿着扫描线$y$从上往下进行扫描。
枚举每一行的间隙，使用上述式子统计贡献，因此我们需要考虑如何快速求出：

$$\sum_{x=x_p+1}^{x_q-1}C_{U[x,y_p]}^kC_{D[x,y_p]}^k$$

可以通过树状数组维护每一个$x$对应的值，前缀和相减即可得到该式。
通过扫描线的移动维护树状数组即可。

---

代码

#include<algorithm>
#include<iostream>
#include<iomanip>
#include<cstring>
#include<cstdlib>
#include<climits>
#include<vector>
#include<cstdio>
#include<cmath>
#include<queue>
using namespace std;

typedef int int_t;
#define int unsigned int

inline const int Get_Int() {
	int num=0,bj=1;
	char x=getchar();
	while(x<'0'||x>'9') {
		if(x=='-')bj=-1;
		x=getchar();
	}
	while(x>='0'&&x<='9') {
		num=num*10+x-'0';
		x=getchar();
	}
	return num*bj;
}

const int maxn=100005,maxk=15;

struct Point {
	int x,y,id;
} p[maxn];

bool cmpx(const Point& a,const Point& b) {
	return a.x<b.x||(a.x==b.x&&a.y<b.y);
}

bool cmpy(const Point& a,const Point& b) {
	return a.y<b.y||(a.y==b.y&&a.x<b.x);
}

int n,k,Left[maxn],Right[maxn],Up[maxn],Down[maxn];

struct BIT {
	int c[maxn];
#define lowbit(x) x&(-x)
	void add(int x,int v) {
		for(int i=x; i<=n; i+=lowbit(i))c[i]+=v;
	}
	int sum(int x) {
		int sum=0;
		for(int i=x; i>=1; i-=lowbit(i))sum+=c[i];
		return sum;
	}
} bit;

int a[maxn],Max=0;
int C[maxn][maxk],ans=0;

int_t main() {
	Get_Int();
	Get_Int();
	n=Get_Int();
	for(int i=1; i<=n; i++) {
		p[i].x=a[i]=Get_Int()+1;
		p[i].y=Get_Int()+1;
		p[i].id=i;
	}
	k=Get_Int();
	sort(a+1,a+n+1);
	int cnt=unique(a+1,a+n+1)-a-1;
	for(int i=1; i<=n; i++)p[i].x=lower_bound(a+1,a+cnt+1,p[i].x)-a;
	sort(p+1,p+n+1,cmpx);
	for(int i=2; i<=n; i++)
		if(p[i].x==p[i-1].x)Down[p[i].id]=Down[p[i-1].id]+1,Max=max(Max,Down[p[i].id]);
	for(int i=n-1; i>=1; i--)
		if(p[i].x==p[i+1].x)Up[p[i].id]=Up[p[i+1].id]+1,Max=max(Max,Up[p[i].id]);
	sort(p+1,p+n+1,cmpy);
	for(int i=2; i<=n; i++)
		if(p[i].y==p[i-1].y)Left[p[i].id]=Left[p[i-1].id]+1,Max=max(Max,Left[p[i].id]);
	for(int i=n-1; i>=1; i--)
		if(p[i].y==p[i+1].y)Right[p[i].id]=Right[p[i+1].id]+1,Max=max(Max,Right[p[i].id]);
	for(int i=0; i<=Max+1; i++) {
		C[i][0]=1;
		for(int j=1; j<=min(k,i); j++)C[i][j]=C[i-1][j-1]+C[i-1][j];
	}
	for(int i=1,j; j=i,i<=n; i=j+1) {
		bit.add(p[i].x,C[Down[p[i].id]+1][k]*C[Up[p[i].id]][k]-C[Down[p[i].id]][k]*C[Up[p[i].id]+1][k]);
		while(p[j+1].y==p[j].y) {
			j++;
			ans+=C[Left[p[j-1].id]+1][k]*C[Right[p[j].id]+1][k]*(bit.sum(p[j].x-1)-bit.sum(p[j-1].x));
			bit.add(p[j].x,C[Down[p[j].id]+1][k]*C[Up[p[j].id]][k]-C[Down[p[j].id]][k]*C[Up[p[j].id]+1][k]);
		}
	}
	printf("%d\n",ans&0x7fffffff);
	return 0;
}