题意

• 传送门

\(N \times M\)的矩形，格点是共\(W\)棵常青树或墓地。对于一块墓地，它的虔诚度为让它正上下左右各恰有\(k\)棵常青树的方法数量。求出整个矩形公墓的虔诚度总和。

对于 \(30\%\) 的数据，满足 \(1 ≤ N, M ≤ 10^3\)。

对于 \(60\%\) 的数据，满足 \(1 ≤ N, M ≤ 10^6\)。

对于 \(100\%\) 的数据，满足 \(1 ≤ N, M ≤ 10^9\)，\(0 ≤ x_i ≤ N\)，\(0 ≤ y_i ≤ M\)，\(1 ≤ W ≤ 10^5\)，\(1 ≤ k ≤ 10\)。

存在 \(50\%\) 的数据，满足 \(1 ≤ k ≤ 2\)。

存在 \(25\%\) 的数据，满足 \(1 ≤ W ≤ 10^4\)。

思路

\(W\)和\(k\)的数据范围较小，考虑在它们身上做文章。

对于\(W\)：离散化，此后最多有\(2*W\)个坐标。
对于\(k\)相关的统计：预处理时间不够，考虑边算边统计。

再整理一遍思路：

1. 按y坐标离散化，并统计每一列一共有多少常青树//，并按列扔进树状数组里

2. 按x坐标离散化，并统计每一行一共有多少常青树

3. 当需要计算时，这时扫描常青树的顺序恰好是按行从左到右的顺序

4. 对于相邻两棵常青树，很容易知道它们之间的空地（墓地）都对应着哪一些列，以及这一些列对于这一块墓地的\(C_{up}^{k} \times C_{down}^{k}\)。
   这个地方是不是得优化……
   每次扫描到一棵树，就可以更新下一行这一列的树的个数了。\(lst_{横坐标}++\)后，把原数\(a\)更新成\(a*(x+1)(y-x-k-1)/(x+1-k)(y-x-1)\).lst记录到现在（包括现在）每一列有多少树。

还要优化，因为没法开1e9的数组。
离散化之后的坐标也要用上。因此，如果直接乘会出问题。加上我们知道每个可能成为十字架的列的前置常青树情况，所以可以通过加减来改变。
在lst[列号]++之前：
给这个坐标加上\(c[lst[ lh ]+1][k]*c[tot[列号]-lst[列号]-2]-c[lst[列号][k]*c[tot[列号]-lst[列号]-1]\).
树状数组：单点修改，区间查询。

5. 把相邻的常青树之间的墓地的\(ans+=C_{left}^k*C_{right}^k*\sum\limits_{这一棵树同一行的上一棵树}^{这一棵树之前}{c[up][k]*c[down][k]}\)，其中left可以现场扫，每一行一共有多少常青树已经记录过了

代码

AC 1.43s

//n m w x y k
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int N=1e5+5;
const ll M=2147483648;
ll n,m;
struct evergreentree {
	ll x,y;
#define x(i) et[i].x
#define y(i) et[i].y
} et[N];
int w,k;
ll lst[N],sx[N],sy[N],zh=1,zl=1;//每一行常青树的总数/列,总行数，总列数
ll c[N][15],tr[N],ans=0;

inline int lowbit(int i){ return i&(-i); }

inline void ad(int p,ll chg){
	while(p<N){
		tr[p]=(tr[p]+chg)%M;
		p+=lowbit(p);
	}
}

inline ll qs(int i){
	ll ret=0;
	while(i){
		ret=(ret+tr[i])%M;
		i-=lowbit(i);
	}
	return ret;
}

inline bool cmpx(evergreentree xx, evergreentree yy) {
	if(xx.x==yy.x) return xx.y<yy.y;
	return xx.x<yy.x;
}

inline bool cmpy(evergreentree xx, evergreentree yy) {
	if(xx.y==yy.y) return xx.x<yy.x;
	return xx.y<yy.y;
}

inline void ini(){
	c[0][0]=1;
	for(int i=1;i<N;i++){
		c[i][0]=1;
		for(int j=1;j<15;j++) c[i][j]=(c[i-1][j]+c[i-1][j-1])%M;
	}
}

int main() {
//	freopen("1.in","r",stdin);
//	
	ini();
	scanf("%lld%lld%d",&n,&m,&w);
	for(int i=1; i<=w; i++) scanf("%lld%lld",&et[i].x,&et[i].y);
	scanf("%d",&k);
	
	sort(et+1,et+w+1,cmpy);
	for(int i=2; i<=w; i++) {
		if(y(i-1)!=y(i)) {
			y(i-1)=zl, sy[zl]++, ++zl;
			continue;
		}
		y(i-1)=zl, sy[zl]++;
	}
	y(w)=zl, sy[zl]++;
	sort(et+1,et+w+1,cmpx);
	for(int i=2; i<=w; i++) {
//		for(int i=1;i<=w;i++) printf("%d ",et[i-1].x);
//		printf("\n");
//		printf("%d %d\n",x(i),y(i));
		if(x(i-1)!=x(i)) {
			x(i-1)=zh, sx[zh]++, ++zh;
			continue;
		}
		x(i-1)=zh, sx[zh]++;
	}
	x(w)=zh, sx[zh]++; 
	
	ll p=0;
	for(int i=1; i<=w; i++) {	
		if(i>1 && x(i-1)==x(i)){
			p++;
			ll m1=c[p][k]%M*c[sx[x(i)]-p][k]%M;
			ll m2=(qs(y(i)-1)-qs(y(i-1))+M)%M;//?
			ans+=m1*m2;
			ans%=M;
		} else p=0;
		lst[y(i)]++;
		ll chg=(c[lst[y(i)]][k]*c[sy[y(i)]-lst[y(i)]][k]-c[lst[y(i)]-1][k]*c[sy[y(i)]-lst[y(i)]+1][k])%M;
		ad(y(i),chg);
	}
	printf("%lld",(ans+M)%M);
	return 0;
}
/*
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
*/