题目传送门：https://www.luogu.org/problemnew/show/P2822 ↓题目大意↓
数据的极限范围：n,m≤2000，k≤21，数据组数≤10000。
由于此题k不大于21，故在计算组合数Cij时，并不需要存储它的真实数值，只需要存储其≤19的所有素因子的个数，判断Cij是否为k的倍数，仅需要判断Cij中各素因子的个数是否大于等于k中的个数即可。基于组合数的性质，我们如果要求出Cij,我们可以通过Ci(j-1)乘上i-j+1然后再除以j即可得到。
下面来考虑如何乘以或除以一个数x。若需要在Ci(j-1)的基础上乘以x，可以考虑将x分解质因数，仅将其≤19的全部素因子与Ci(j-1)的素因子个数进行累加。除法同理，加法改成减法即可。
最后维护一个二维数组b。若b[i][j]=1，则表示Cij是k的倍数。输入n,m时，将b[1..n][1..m]进行累加即可得出答案。很明显这么操作依然会TLE，使用另一数组维护b[i][j]的前缀和即可。
时间复杂度为O(n*m+T)。 但常数很大（本地均为0.3s左右）。
#include<iostream>
#include<cstdio>
#include<cstring>
#define M 2000
using namespace std;
int p[]={,,,,,,,};
struct cg{
int a[];
cg(){memset(a,,sizeof(a));}
cg(int x){
for(int i=;i<;i++)
while(x%p[i]==) a[i]++,x/=p[i];
}
friend cg operator *(cg a,int x){
for(int i=;i<;i++)
while(x%p[i]==) a.a[i]++,x/=p[i];
return a;
}
friend cg operator /(cg a,int x){
for(int i=;i<;i++)
while(x%p[i]==) a.a[i]--,x/=p[i];
return a;
}
friend bool operator +(cg a,int x){
cg c=a;
for(int i=;i<;i++)
while(x%p[i]==){
c.a[i]--;x/=p[i];
if(c.a[i]<) return ;
}
return ;
}
}a[M+][M+];
int b[M+][M+]={},k; void init(){
for(int i=;i<=M;i++){
int zhi=i>>;
for(int j=;j<=zhi;j++)
a[i][j]=a[i][i-j]=a[i][j-]*(i-j+)/j;
}
for(int i=;i<=M;i++){
int zhi=i>>;
for(int j=;j<=zhi;j++){
b[i][j]=b[i][i-j]=a[i][j]+k;
}
}
for(int i=;i<=M;i++)
for(int j=;j<=M;j++){
b[i][j]=b[i-][j]+b[i][j-]-b[i-][j-]+b[i][j];
}
} int main(){
freopen("problem.in","r",stdin);
freopen("problem.out","w",stdout);
int cas; cin>>cas>>k;
init();
while(cas--){
int x,y; scanf("%d%d",&x,&y);
printf("%d\n",b[x][y]);
}
}