BZOJ3434: [Wc2014]时空穿梭

我们考虑怎么选取c个点满足有直线同时经过他们
设$A_i$表示第i个点，用向量$b_i=(x1,x2,x3...xn)$表示$A_i+1$与$A_i$的差
由于每一维的$vi（xi=pi+tvi）$是不变的，对于相邻的两个点，$\Delta t$是相同的，所以任一向量$bi=(x1,x2,x3...xn)$，$x_i$之间的比值是不变的，即$v_i$之间的比值，所以对于任一向量$bi=(x1,x2,x3...xn)$，都可以将它写成$bi=(\Delta t*v1,\Delta t*v2...\Delta t*vn)$，将$v_i$之间的比值写成整数$a_i$之间的比值，得$bi=(k*c1,k*c2,...k*cn)$

我们考虑枚举每一维的极差$a1,a2,...an$，设$d=gcd(a1,a2,...an)$，则$a_i=c_i*d$，则将$d$分给这$c-1$个$bi$，共有$C_{d-1}^{c-2}$种分法，原理是插板法，即因为每一维严格递增，所以每个$b_i$分到的不能为0，相当于$1\backsim d$，之间有$d-1$个空位，在这些中选出$c-2$个位置插隔板代表这$c-1$个数之间的分界线（感觉没讲清楚，意会一下？）
因为所有坐标为正整数，第$i$维坐标不超过$M_i$，所以极差不超过$M_i-1$
设$m_i=M_i-1$
所以有

$$ans=\sum_{a_1=1}^{m_1}(M_1-a_1)\sum_{a_2=1}^{m_2}(M_2-a_2)...\sum_{a_n=1}^{m_n}(M_n-a_n)C_{gcd(a1,a2...an)-1}^{c-2}$$
$=\sum_{d}C_{d-1}^{c-2}\prod_{i=1}^n\sum_{a_i=1}^{\lfloor \frac{m_i}{d}\rfloor}(M_i-a_id)\ \epsilon ((a1,a2,...an)=1)$
反演一下

$=\sum_{d}C_{d-1}^{c-2}\prod_{i=1}^n\sum_{a_i=1}^{\lfloor \frac{m_i}{d}\rfloor}(M_i-a_id)\sum_{D|(a1,a2...an)}μ(D)$

$=\sum_{d}C_{d-1}^{c-2}\sum_{D}μ(D)\prod_{i=1}^n\sum_{a_i=1}^{\lfloor \frac{m_i}{Dd}\rfloor}(M_i-a_iDd)$

枚举$D=Dd$
$=\sum_{D}\sum_{d|D}C_{d-1}^{c-2}μ(\frac{D}{d})\prod_{i=1}^n\sum_{a_i=1}^{\lfloor \frac{m_i}{D}\rfloor}(M_i-a_iD)$

可以发现$f(D,c)=\sum_{d|D}C_{d-1}^{c-2}μ(\frac{D}{d})$这个东西只与$D$与$c$有关，$c$的范围不大，可以对所有的$D$和$c$预处理出$f(D,c)$（枚举$c$后，对每个$d$，枚举它的倍数$D$，根据调和级数，这样总的复杂度是$O(cmlogm)$）

到这里，因为后面的柿子和D仍然相关，仍然需要枚举$m$个$D$，因为$m$很大，而且有多组数据，观察到使每个$\lfloor \frac{m_i}{D}\rfloor$相同的$D$最多有$O(n\sqrt m)$段，这个复杂度可以接受，所以考虑继续化，

观察$\prod_{i=1}^n\sum_{a_i=1}^{\lfloor \frac{m_i}{D}\rfloor}(M_i-a_iD)$这个东西，可以化成
$\prod_{i=1}^nMi\lfloor \frac{m_i}{D}\rfloor-\dfrac{(D+D\lfloor \frac{m_i}{D}\rfloor)\lfloor \frac{m_i}{D}\rfloor}{2}$

$=\prod_{i=1}^nMi\lfloor \frac{m_i}{D}\rfloor-D\dfrac{(1+\lfloor \frac{m_i}{D}\rfloor)\lfloor \frac{m_i}{D}\rfloor}{2}$

接着把这$n$个多项式相乘$O(n^2)$展开成形如$a_0+a_1D+a_2D^2+.....a_nD^n$的多项式，对于$D^i$的系数$a_i$，只与$\lfloor \frac{m_i}{D}\rfloor$有关

所以在维护了$g(D,c,j)=f(D,c)D^j$的前缀和后，可以将原式按$\lfloor \frac{m_i}{D}\rfloor$分成$O(n\sqrt m)$段，每段内暴力$O(n^2)$展开多项式各项系数，乘上前面的柿子这一段的和

时间复杂度$O(nmc+cmlogm+Tn^3\sqrt m)$

code：

#include<set>
#include<map>
#include<deque>
#include<queue>
#include<stack>
#include<cmath>
#include<ctime>
#include<bitset>
#include<string>
#include<vector>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<climits>
#include<complex>
#include<iostream>
#include<algorithm>
#define ll long long
using namespace std;

inline void down(int &x,const int y){if(x>y)x=y;}
const int maxn = 12;
const int maxc = 21;
const int maxm = 110000;
const int Mod = 10007;

int p[maxm],pri,miu[maxm];
bool v[maxm];
int C[maxm][maxc],f[maxc][maxm],pg[maxc][maxm][maxn];
//                 miu*C            pre G(i,c)×i^j
void pre()
{
    C[0][0]=1;
    for(int i=1;i<maxm;i++)
    {
        C[i][0]=1;
        for(int j=1;j<maxc;j++) C[i][j]=(C[i-1][j-1]+C[i-1][j])%Mod;
    }

    miu[1]=1; ///
    for(int i=2;i<maxm;i++)
    {
        if(!v[i]) miu[i]=-1,p[++pri]=i;
        for(int j=1;j<=pri;j++)
        {
            int k=i*p[j];
            if(k>=maxm) break;
            v[k]=true;
            if(i%p[j]==0) break;
            miu[k]=-miu[i];
        }
    }

    for(int c=2;c<maxc;c++)
    {
        for(int i=1;i<maxm;i++)
        {
            for(int j=i,ji=1;j<maxm;j+=i,ji++)
                (f[c][j]+=C[i-1][c-2]*miu[ji]%Mod)%=Mod;
        }
    }

    for(int c=2;c<maxc;c++)
    {
        for(int i=1;i<maxm;i++)
        {
            for(int j=0,ij=1;j<maxn;j++,ij=ij*i%Mod)
            {
                pg[c][i][j]=(pg[c][i-1][j]+f[c][i]*ij%Mod)%Mod;
            }
        }
    }
}

int n,c;
int m[maxn];
int a[maxn],an;

int main()
{
    pre();

    int t; scanf("%d",&t);
    while(t--)
    {
        int mn=maxm;
        scanf("%d%d",&n,&c);
        for(int i=1;i<=n;i++) scanf("%d",&m[i]),m[i]--,down(mn,m[i]);

        int ret=0;
        for(int i=1;i<=mn;)
        {
            int r=mn;
            for(int j=1;j<=n;j++) down(r,m[j]/(m[j]/i));

            a[an=0]=1;
            for(int j=1;j<=n;j++)
            {
                int mi=(m[j]/i)%Mod;
                int x1=(m[j]+1)*mi%Mod,x2=-mi*(mi+1)/2%Mod;
                a[an+1]=a[an]*x2%Mod;
                for(int l=an;l>=0;l--) a[l]=(a[l]*x1%Mod+(l?a[l-1]:0)*x2%Mod)%Mod;
                an++;
            }
            for(int j=0;j<=an;j++) (ret+=a[j]*((pg[c][r][j]-pg[c][i-1][j])%Mod)%Mod)%=Mod;

            i=r+1;
        }
        if(ret<0) ret+=Mod;
        printf("%d\n",ret);
    }

    return 0;
}