【WC2014】时空穿梭【组合数】【莫比乌斯反演】【整除分块】【暴力多项式】

题意：$T$ 组数据，给一个 $n$ 维空间，第 $i$ 维大小为 $[1,m_i]\cap \mathbb{Z}$，求大小为 $c$ 的严格偏序上升的共线点集个数。答案模 $10007$。

$T\le 100,n\le 11,m\le 10^5,c\le 20$

神仙题

设空间大小向量为 $m$,$M$ 为坐标最大值。下文中若未特殊说明，对向量无公开定义的运算均定义为其每一维坐标的运算。

考虑第一个点到最后一个点的差向量 $x$，它中间有 $\gcd (x)-1$ 个空隙可以放，$\gcd$ 为所有坐标数值的 $\gcd$。所以

$$ans=\sum _{x\le m}\left(\genfrac{}{}{0ex}{}{\gcd (x)-1}{c-2}\right)\prod _{i=1}^n(m_i-x_i)$$

看到有个 $\gcd$ ，想到反演。

$$\begin{gathered} ans=\sum _{d=1}^M\left(\genfrac{}{}{0ex}{}{d-1}{c-2}\right)\sum _{x\le m}[\gcd (x)=d]\prod _{i=1}^n(m_i-x_i) \\ =\sum _{d=1}^M\left(\genfrac{}{}{0ex}{}{d-1}{c-2}\right)\sum _{x\le \frac{m}{d}}[\gcd (x)=1]\prod _{i=1}^n(m_i-d{x}_i) \\ =\sum _{d=1}^M\left(\genfrac{}{}{0ex}{}{d-1}{c-2}\right)\sum _{x\le \frac{m}{d}}\sum _{k\mid x}\mu (k)\prod _{i=1}^n(m_i-d{x}_i) \\ =\sum _{d=1}^M\left(\genfrac{}{}{0ex}{}{d-1}{c-2}\right)\sum _{k\mid x}\mu (k)\sum _{x\le \frac{m}{dk}}\prod _{i=1}^n(m_i-dk{x}_i) \\ =\sum _{d=1}^M\sum _{x\le \frac{m}{d}}\prod _{i=1}^n(m_i-d{x}_i)\sum _{k\mid d}\left(\genfrac{}{}{0ex}{}{k-1}{c-2}\right)\mu (\frac{d}{k}) \end{gathered}$$

令 $g(n)=\sum _{d\mid n}\left(\genfrac{}{}{0ex}{}{d-1}{c-2}\right)\mu (\frac{n}{d})$,这个随便预处理一下就能算出来。

$$\begin{gathered} ans=\sum _{d=1}^M\sum _{x\le \frac{m}{d}}\prod _{i=1}^n(m_i-d{x}_i)g(d) \\ =\sum _{d=1}^{M}g(d)\prod _{i=1}^n\sum _{x=1}^{\lfloor \frac{m_i}{d}\rfloor }(m_i-dx) \end{gathered}$$

这样大概是 $\mathrm{O}(Tmn\log n)$ 的，有点卡。

看到整除，立刻想到整除分块。

注意到后面是一个关于 $d$ 的 $n$ 次多项式（每个和式次数为 $1$），如果要整除分块的话只需要求这个多项式的前缀和。

设这个多项式系数为 $a_i$,那么一次区间查询

$$\begin{gathered} \sum _{d=l}^{r}g(d)f(d) \\ =\sum _{d=l}^{r}g(d)\sum _{i=0}^n{a}_i{d}^i \\ =\sum _{i=0}^n{a}_i\sum _{d=l}^{r}g(d)d^i \end{gathered}$$

后面是和询问无关的，直接预处理就可以 $\mathrm{O}(n)$ 回答。

现在的问题是如何维护多项式。 整除分块时一共变化了 $\mathrm{O}(n\sqrt{m})$ 次，每次变化的时候把操作的 $1$ 次多项式做一个多项式乘法或除法就可以了。

然后是喜闻乐见的除 $0$ 问题。记下除了几个 $0$ 多项式即可，但其他人都没提到这个，代码也看不懂，一脸懵逼，可能写法不太一样吧……

复杂度大概是 $\mathrm{O}(ncm+T{n}^2\sqrt{m})$

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cctype>
using namespace std;
typedef long long ll;
inline int read()
{
	int ans=0;
	char c=getchar();
	while (!isdigit(c)) c=getchar();
	while (isdigit(c)) ans=(ans<<3)+(ans<<1)+(c^48),c=getchar();
	return ans;
}
const int MOD=10007;
inline int add(const int& x,const int& y){return x+y>=MOD? x+y-MOD:x+y;}
inline int dec(const int& x,const int& y){return x<y? x-y+MOD:x-y;}
inline int qpow(int a,int p)
{
	int ans=1;
	a%=MOD;
	while (p)
	{
		if (p&1) ans=ans*a%MOD;
		a=a*a%MOD,p>>=1;
	}
	return ans;
}
const int N=1e5,MAXN=N+5;
int np[MAXN],pl[MAXN],cnt;
inline void init()
{
	np[1]=1;
	for (int i=2;i<=N;i++)
	{
		if (!np[i]) pl[++cnt]=i;
		for (int j=1,x;(x=i*pl[j])<=N;j++)
		{
			np[x]=1;
			if (i%pl[j]==0) break;
		}
	}
}
int C[MAXN][20],g[22][MAXN],S[12][22][MAXN];
int a[12],n,c,lim[12],m[12],zero;
inline void mul(int x,int y)
{
	if (!x&&!y) return (void)(++zero);
	for (int i=n;i>=1;i--) a[i]=(a[i-1]*x+a[i]*y)%MOD;
	a[0]=a[0]*y%MOD;
}
inline void div(int x,int y)
{
	if (!y)
	{
		if (!x) return (void)(--zero);
		int ix=qpow(x,MOD-2);
		for (int i=0;i<n;i++) a[i]=a[i+1]*ix%MOD;
		a[n]=0; 
		return;
	}
	int iy=qpow(y,MOD-2);
	a[0]=a[0]*iy%MOD;
	for (int i=1;i<=n;i++) a[i]=dec(a[i],a[i-1]*x%MOD)*iy%MOD;
}
inline int calc(int d,int l,int r)
{
	for (int i=1;i<=n;i++)
	{
		int cur=m[i]/d;
		if (cur<lim[i])
		{
			div(dec(0,lim[i]*(lim[i]+1ll)/2%MOD),(ll)lim[i]*m[i]%MOD);
			lim[i]=cur;
			mul(dec(0,lim[i]*(lim[i]+1ll)/2%MOD),(ll)lim[i]*m[i]%MOD);
		}
	}
	if (zero) return 0;
	int ans=0;
	for (int i=0;i<=n;i++) 
		ans=(ans+a[i]*dec(S[i][c][r],S[i][c][l-1]))%MOD;
	return ans;
}
int main()
{
	init();
	C[0][0]=1;
	for (int i=1;i<=N;i++)
	{
		C[i][0]=1;
		for (int j=1;j<20;j++) C[i][j]=(C[i-1][j]+C[i-1][j-1])%MOD;
	}
	for (int t=2;t<22;t++)
	{
		for (int i=1;i<=N;i++) g[t][i]=C[i-1][t-2];
		for (int j=1;j<=cnt;j++)
			for (int i=N/pl[j];i>=1;i--)
				g[t][i*pl[j]]=dec(g[t][i*pl[j]],g[t][i]);
		for (int i=0;i<12;i++)
			for (int d=1;d<=N;d++)
				S[i][t][d]=add(S[i][t][d-1],qpow(d,i)*g[t][d]%MOD);			
	}
	for (int T=read();T;T--)
	{
		zero=0;
		n=read(),c=read();
		int M=N;
		for (int i=1;i<=n;i++) M=min(M,lim[i]=m[i]=read());
		memset(a,0,sizeof(a));
		a[0]=1;
		for (int i=1;i<=n;i++) mul(dec(0,m[i]*(m[i]+1ll)/2%MOD),(ll)m[i]*m[i]%MOD);
		int ans=0;
		for (int l=1,r;l<=M;l=r+1)
		{
			r=M;
			for (int i=1;i<=n;i++) r=min(r,m[i]/(m[i]/l));
			ans=add(ans,calc(l,l,r));
		}
		printf("%d\n",ans);
	}
	return 0;
}