bzoj3434 [Wc2014]时空穿梭（莫比乌斯函数+多项式乘法）

蒟蒻我第一步共线就没想出来x暴力dp T飞。直接爆零qaq
我们记C个点分别为$C_i$,V为$\{v_1,v_2...v_n\}$,令$B_i=C_i-C_{i-1}$。则我们发现$B_i=(t_i-t_{i-1})V$，也就是说每一个B的各维坐标都是成比例的。我们设第i维的坐标极差为$\Delta x_i$（$C_n-C_1$）,，则这条直线上最多有（除了起终点）$gcd(\Delta x_i)-1$个点。我们要从中再选出c-2个点。而这样的直线起点有$(M_i-\Delta x_i)$种。显然每一维的极差最大为$M_i-1$，我们设为$m_i$，并令$m=min\{m_i\}$所以我们枚举所有的极差，得到

$Ans=\sum\limits_{\Delta x_1=1}^{m_1}\sum\limits_{\Delta x_2=1}^{m_2}...\sum\limits_{\Delta x_n=1}^{m_n}C_{gcd(x_1,x_2,...,x_n)-1}^{c-2}\prod\limits_{i=1}^n(M_i-\Delta x_i)$

设$d=gcd(x_1,x_2,...,x_n)$,则

$Ans=\sum\limits_{d=1}^mC_{d-1}^{c-2}\sum\limits_{\Delta x_1=1}^{\lfloor\frac{m_1}{d}\rfloor}\sum\limits_{\Delta x_2=1}^{\lfloor\frac{m_2}{d}\rfloor}...\sum\limits_{\Delta x_n=1}^{\lfloor\frac{m_n}{d}\rfloor}[gcd(x_1,x_2,...,x_n)==1]\prod\limits_{i=1}^n(M_i-\Delta x_i)$

莫比乌斯反演，得到

$Ans=\sum\limits_{d=1}^mC_{d-1}^{c-2}\sum\limits_{\Delta x_1=1}^{\lfloor\frac{m_1}{d}\rfloor}\sum\limits_{\Delta x_2=1}^{\lfloor\frac{m_2}{d}\rfloor}...\sum\limits_{\Delta x_n=1}^{\lfloor\frac{m_n}{d}\rfloor}\sum\limits_{t|gcd(x_1,...x_n)}\mu(t)\prod\limits_{i=1}^n(M_i-\Delta x_id)$

$Ans=\sum\limits_{d=1}^mC_{d-1}^{c-2}\sum\limits_{t=1}^{\lfloor m/d\rfloor}\mu(t)\sum\limits_{\Delta x_1=1}^{\lfloor\frac{m_1}{dt}\rfloor}\sum\limits_{\Delta x_2=1}^{\lfloor\frac{m_2}{dt}\rfloor}...\sum\limits_{\Delta x_n=1}^{\lfloor\frac{m_n}{dt}\rfloor}\prod\limits_{i=1}^n(M_i-\Delta x_idt)$

令h=dt,
$Ans=\sum\limits_{h=1}^m\sum\limits_{\Delta x_1=1}^{\lfloor\frac{m_1}{h}\rfloor}\sum\limits_{\Delta x_2=1}^{\lfloor\frac{m_2}{h}\rfloor}...\sum\limits_{\Delta x_n=1}^{\lfloor\frac{m_n}{h}\rfloor}\prod\limits_{i=1}^n(M_i-\Delta x_ih)\sum\limits_{d|h}C_{d-1}^{c-2}\mu(h/d)$

$Ans=\sum\limits_{h=1}^m\prod\limits_{i=1}^n\sum\limits_{x=1}^{\lfloor\frac{m_i}{h}\rfloor}(M_i-xh)\sum\limits_{d|h}C_{d-1}^{c-2}\mu(h/d)$

我们令$G[h]=\sum\limits_{d|h}C_{d-1}^{c-2}\mu(h/d)$
显然我们是可以$O(mlogm)$的预处理出G数组的。
令$f(m,h)=\sum\limits_{x=1}^{\lfloor\frac{m}{h}\rfloor}(M_i-xh)$
$f(m,h)=\lfloor\frac{m}{h}\rfloor*M_i-h*(1+\lfloor\frac{m}{h}\rfloor)*\lfloor\frac{m}{h}\rfloor/2$可以O(1)算。

所以
$Ans=\sum\limits_{h=1}^m\prod\limits_{i=1}^nf(m_i,h)G[h]$
复杂度$O(nm)$

我们考虑进一步优化，由于n元组{$m_1/h...m_n/h$}的取值最多$n\sqrt{m}$种，而确定了他们的取值后f函数就是一个形如px+q的多项式，假设在l…r上这n个多项式是相同的，我们可以$O(n^2)$的算出这n个多项式的乘积，记作P(h),则
$Ans_{l...r}=\sum\limits_{h=l}^rP(h)G[h]$
我们记$P(h)=\sum\limits_{k=0}^na_kh^k$则
$Ans_{l...r}=\sum\limits_{k=0}^na_k\sum\limits_{h=l}^rG[h]*h^k$。
我们可以预处理出$G[h]h^k$的前缀和。

所以总的复杂度就是$O(nmc+cmlogm+T*n^3\sqrt{m})$

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <queue>
using namespace std;
#define ll long long
#define inf 0x3f3f3f3f
#define N 100010
#define mod 10007
inline int read(){
    int x=0,f=1;char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-') f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9') x=x*10+ch-'0',ch=getchar();
    return x*f;
}
int tst,n,m[13],G[N],mu[N],C[N][20],sum[22][13][N];
bool notprime[N];int a[13],prime[N],tot=0;
inline void getmu(){
    notprime[1]=1;mu[1]=1;
    for(int i=2;i<=100000;++i){
        if(!notprime[i]){
            prime[++tot]=i;mu[i]=-1;
        }for(int j=1;prime[j]*i<=100000;++j){
            notprime[prime[j]*i]=1;
            if(i%prime[j]==0){
                mu[prime[j]*i]=0;break;
            }mu[prime[j]*i]=-mu[i];
        }
    }for(int i=0;i<=100000;++i) C[i][0]=1;C[1][1]=1;
    for(int i=2;i<=100000;++i)
        for(int j=1;j<=18;++j) C[i][j]=(C[i-1][j-1]+C[i-1][j])%mod;
    for(int c=2;c<=20;++c){
        memset(G,0,sizeof(G));
        for(int i=1;i<=100000;++i)
            for(int j=i;j<=100000;j+=i)
                G[j]=(G[j]+C[i-1][c-2]*mu[j/i])%mod;
        for(int i=1;i<=100000;++i){
            int tmp=1,base=i%mod;
            for(int k=0;k<=11;++k){
                sum[c][k][i]=(sum[c][k][i-1]+G[i]*tmp)%mod;
                tmp=tmp*base%mod;
            }
        }
    }
}
inline void calc(int h){//计算n个多项式(px+q)的乘积
    memset(a,0,sizeof(a));a[0]=1;
    for(int i=1;i<=n;++i){
        int t=m[i]/h;
        int p=(-(ll)(t+1)*t/2)%mod;
        int q=(ll)t*(m[i]+1)%mod;
        for(int j=n;j>=1;--j)
            a[j]=(a[j-1]*p+a[j]*q)%mod;
        a[0]=(ll)a[0]*q%mod;
    }
}
int main(){
//  freopen("space1.in","r",stdin);
//  freopen("space1.out","w",stdout);
    tst=read();getmu();
    while(tst--){
        n=read();int c=read();int mn=inf,last,ans=0;
        for(int i=1;i<=n;++i) m[i]=read()-1,mn=min(mn,m[i]);
        for(int i=1;i<=mn;i=last+1){
            last=inf;for(int j=1;j<=n;++j) last=min(last,m[j]/(m[j]/i));
            calc(i);for(int j=0;j<=n;++j)
            ans=(ans+(ll)a[j]*(sum[c][j][last]-sum[c][j][i-1]))%mod;
        }if(ans<0) ans+=mod;printf("%d\n",ans);
    }return 0;
}