题目描述
小X驾驶着他的飞船准备穿梭过一个n维空间,这个空间里每个点的坐标可以用n个实数来表示,即(x1,x2,...,xn) 。
为了穿过这个空间,小 X 需要在这个空间中选取 c(c≥2) 个点作为飞船停留的地方,而这些点需要满足以下三个条件:
每个点的每一维坐标均为正整数,且第 ii 维坐标不超过 \(m_i\) 。
第 \(i + 1\) (\(1 \leq i < c\))个点的第 j (\(1 \leq j \leq n\)) 维坐标必须严格大于第 i 个点的第 j 维坐标。
存在一条直线经过所选的所有点。在这个 n 维空间里,一条直线可以用 2n个实数 \(p_1, p_2, … , p_n\), \(v_1, v_2, … , v_n\)表示。 直线经过点 (x1,x2,...,xn) ,当且仅当存在实数 t,使得对 \(i =1 … n\) 均满足 \(x_i = p_i + tv_i\)。
小 X 还没有确定他的最终方案,请你帮他计算一下一共有多少种不同的方案满足他的要求。由于答案可能会很大,你只需要输出答案 mod 10 007后的值。
Input
输入文件 space.in 的第一行包含一个正整数 TT ,表示有 TT 组数据需要求解。
每组数据包含两行,第一行包含两个正整数 n, c (\(c \geq 2\)),分别表示空间的维数和需要选择的暂停点个数。
第二行包含 n 个正整数,依次表示 \(m_1\), \(m_2\), … , \(m_n\)。
Output
输出文件 space.out 包含 T 行,每行包含一个非负整数,依次对应每组数据的答案。
Solution
神奇反演题。。
先考虑二维的情况:
考虑枚举左下角和右上角两个点的坐标差值,中间的整点可以随意填,这个矩形还可以在平面上平移,所以答案可以形式化的写成:
\[ ans=\sum_{i=1}^{n}\sum_{j=1}^{m}(n-i)(m-j)\binom{\gcd(i,j)-1}{c-2} \]
然后对\(\gcd(i,j)\)莫比乌斯反演可得:
\[ \begin{align} ans=&\sum_{T=1}^{\min(n,m)}(\lfloor\frac{n}{T}\rfloor\cdot n-T\cdot \frac{\lfloor\frac{n}{T}\rfloor(\lfloor\frac{n}{T}\rfloor+1)}{2})\\& \cdot(\lfloor\frac{m}{T}\rfloor\cdot m-T\cdot \frac{\lfloor\frac{m}{T}\rfloor(\lfloor\frac{m}{T}\rfloor+1)}{2})\sum_{d|T}\mu(\frac{T}{d})\binom{d-1}{c-2} \end{align} \](太长了换个行)
然后发现无论多少维最后一个求和都是一样的,所以可以推广到\(n\)维,设每一维的限制为\(m_i\),最小的为\(mn\),可得:
\[ ans=\sum_{T=1}^{mn}\prod_{i=1}^n(\lfloor\frac{m_i}{T}\rfloor\cdot m_i-T\cdot \frac{\lfloor\frac{m_i}{T}\rfloor(\lfloor\frac{m_i}{T}\rfloor+1)}{2})\sum_{d|T}\mu(\frac{T}{d})\binom{d-1}{c-2} \](然而这并不能过掉此题。。)
这样复杂度是\(O(Tnm)\)左右。
然后考虑下中间那个连乘项,发现中间有向下取整,且\(n\)很小,可以尝试丧心病狂的进行数论分块。
假设当前数论分块的区间所有的\(\lfloor\frac{m_i}{T}\rfloor\)都不变,那么可以把那一项看作是很多个关于\(T\)的一次多项式相乘,其他的都是常量,然后暴力展开,设系数为\(a_i\),中间可以写成:
\[ \sum_{i=0}^{n}a_iT^i \]
设后面一项为\(g(i)\),即:
\[ g(n)=\sum_{d|n}\mu(\frac{n}{d})\binom{d-1}{c-2} \]
答案可以化成:
\[ ans=\sum_{T=1}^{mn}\sum_{i=0}^{n}a_iT^ig(T)=\sum_{i=0}^{n}a_i\sum_{T=1}^{mn}T^ig(T) \]
然后预处理后面那一项,系数暴力算,就行了。
注意\(c\)只有20种取值,可以一开始预处理出\(g(c,T)\),复杂度\(O(cm\log m)\),总复杂度\(O(cm\log m+n^2m+Tn^2\sqrt{m})\)
#include<bits/stdc++.h>
using namespace std;
#ifdef ONLINE_JUDGE
#define getchar() ((p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin)),p1==p2)?EOF:*p1++)
#endif
namespace fast_IO {
char buf[1<<21],*p1=buf,*p2=buf;
template <typename T> inline void read(T &x) {
x=0;T f=1;char ch=getchar();
for(;!isdigit(ch);ch=getchar()) if(ch=='-') f=-f;
for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';x*=f;
}
template <typename T,typename... Args> inline void read(T& x,Args& ...args) {
read(x),read(args...);
}
char buf2[1<<21],a[80];int p,p3=-1;
inline void flush() {fwrite(buf2,1,p3+1,stdout),p3=-1;}
template <typename T> inline void write(T x) {
if(p3>(1<<20)) flush();
if(x<0) buf2[++p3]='-',x=-x;
do {a[++p]=x%10+48;} while(x/=10);
do {buf2[++p3]=a[p];} while(--p);
buf2[++p3]='\n';
}
template <typename T,typename... Args> inline void write(T x,Args ...args) {
write(x),write(args...);
}
}
using fast_IO :: read;
using fast_IO :: write;
using fast_IO :: flush;
const int mod = 10007;
const int inf = 1e9;
const int maxn = 1e5+10;
int pri[maxn],vis[maxn],g[21][maxn],mu[maxn],c[maxn][21],tot,sum[12][21][maxn],pw[maxn][13];
int qpow(int a,int x) {
int res=1;
for(;x;x>>=1,a=1ll*a*a%mod) if(x&1) res=res*a%mod;
return res;
}
void prepare() {
mu[1]=1;
for(int i=2;i<maxn;i++) {
if(!vis[i]) pri[++tot]=i,mu[i]=-1;
for(int j=1;j<=tot&&i*pri[j]<maxn;j++) {
vis[i*pri[j]]=1;
if(i%pri[j]==0) break;
mu[i*pri[j]]=-mu[i];
}
}
c[0][0]=1;
for(int i=1;i<maxn;i++) {
c[i][0]=1;
for(int j=1;j<=20;j++) c[i][j]=(c[i-1][j-1]+c[i-1][j])%mod;
}
for(int i=2;i<=20;i++)
for(int d=1;d<maxn;d++)
for(int T=d;T<maxn;T+=d)
(g[i][T]+=mu[T/d]*c[d-1][i-2])%=mod;
for(int j=1;j<=20;j++)
for(int k=1;k<maxn;k++)
for(int i=0,p=1;i<=11;i++,p=p*k%mod)
sum[i][j][k]=(sum[i][j][k-1]+g[j][k]*p%mod)%mod;
}
int m[20],C,a[12],b[12],n;
void calc(int T) {
memset(a,0,sizeof a);
int t;a[0]=1;
for(int i=1;i<=n;i++) {
t=m[i]/T;
int x=1ll*t*m[i]%mod,y=-1ll*t*(t+1)/2%mod;
for(int j=i;j;j--) a[j]=(a[j]*x+a[j-1]*y)%mod;
a[0]*=x;a[0]%=mod;
}
}
void solve() {
read(n,C);for(int i=1;i<=n;i++) read(m[i]);
int mn=inf;for(int i=1;i<=n;i++) mn=min(mn,m[i]);
int T=1,ans=0;
while(T<=mn) {
int pre=T;T=inf;
for(int i=1;i<=n;i++) T=min(m[i]/(m[i]/pre),T);
calc(T);
for(int i=0;i<=n;i++)
ans=(ans+a[i]*(sum[i][C][T]-sum[i][C][pre-1])%mod)%mod;
T++;
}
write((ans%mod+mod)%mod);
}
int main() {
prepare();//cerr << (double)clock()/1e6 << endl;
int t;read(t);
while(t--) solve();
flush();
return 0;
}