[数学记录]CF961G Partitions

command_block · 2019-07-20 22:58:56 · 个人记录

CF961G Partitions

讲真这个翻译我没看懂……

给出n个物品,每个物品有一个权值w_i

定义一个集合S的权值W(S)=|S|\sum\limits_{x\in S}w_x

把全集划分为k个集合S_{1...k},那么贡献是\sum\limits_{i=1}^kW(S_i)

求所有划分方法的贡献和\ \bmod\ 10^9+7.

n,k\le2\times10^5,w_i\le10^9

首先,能想到:所有物品在划分面前都是平等的。

那么每个物品的贡献常数都是一样的,设为p,那么答案是p\sum\limits_{i=1}^nw_i

我们考虑枚举包含该物品的集合大小:

p=\sum\limits_{i=1}^ni\dbinom{n-1}{i-1}\begin{Bmatrix}n-i\\k-1\end{Bmatrix}

(大括号括起来的是第二类斯特林数,关于斯特林数)

乘的那个i是因为W(S)有系数|S|.

这东西就是快速求一列第二类斯特林数?可是模数非NTT模数,咕咕咕?

一个方法是MTT卡常,之前没什么脑子只会拖板子,所以有如下代码 :

#include<algorithm>
#include<cstdio>
#include<cmath>
#define mod 1000000007
#define lim 32000
#define Maxn 200500
#define ll long long
using namespace std;
const long double Pi=acos(-1);
int powM(long long a,int t=mod-2)
{
  long long ans=1,buf=a;
  while(t){
    if(t&1)ans=(ans*buf)%mod;
    buf=(buf*buf)%mod;
    t>>=1;
  }return ans;
}
int n,k;
ll fac[Maxn],inv[Maxn];
int tr[Maxn<<2];
struct CP
{
  long double x,y;
  CP operator + (const CP& B) const
  {return (CP){x+B.x,y+B.y};}
  CP operator - (const CP& B) const
  {return (CP){x-B.x,y-B.y};}
  CP operator * (const CP& B) const
  {return (CP){x*B.x-y*B.y,x*B.y+y*B.x};}
}P1[Maxn<<2],P2[Maxn<<2],Q[Maxn<<2];
void FFT(CP *f,int op,int n)
{
  for (int i=0;i<n;i++)
   if (i<tr[i])swap(f[i],f[tr[i]]);
  for(int p=2;p<=n;p<<=1){
    int len=p>>1;
    CP tmp=(CP){std::cos(Pi/len),op*std::sin(Pi/len)};
    for(int k=0;k<n;k+=p){
      CP buf=(CP){1,0};
      for(int l=k;l<k+len;l++){
        CP tt=buf*f[len+l];
        f[len+l]=f[l]-tt;
        f[l]=f[l]+tt;
        buf=buf*tmp;
      }//F(x)=FL(x^2)+x*FR(x^2)
       //F(W^k)=FL(w^k)+W^k*FR(w^k)
       //F(W^{k+n/2})=FL(w^k)-W^k*FR(w^k)
    }
  }
}
void times(ll *f,ll *g,int len,int siz)
{
  for (int i=0;i<len;i++){
    P1[i]=(CP){f[i]/lim,f[i]%lim};
    P2[i]=(CP){f[i]/lim,-f[i]%lim};
    Q[i]=(CP){g[i]/lim,g[i]%lim};
  }int sav=len+len;len=1;
  for (;len<sav;len<<=1);
  for (int i=1;i<len;i++)
    tr[i]=tr[i>>1]>>1|((i&1)?len>>1:0);
  FFT(P1,1,len);FFT(P2,1,len);FFT(Q,1,len);
  for (int i=0;i<len;i++){Q[i].x/=len;Q[i].y/=len;}
  for (int i=0;i<len;i++)P1[i]=P1[i]*Q[i];
  for (int i=0;i<len;i++)P2[i]=P2[i]*Q[i];
  FFT(P1,-1,len);FFT(P2,-1,len);
  for (int i=0;i<siz;i++){
    ll a1b1,a1b2,a2b1,a2b2;
    a1b1=(ll)floor((P1[i].x+P2[i].x)/2+0.49)%mod;
    a1b2=(ll)floor((P1[i].y+P2[i].y)/2+0.49)%mod;
    a2b1=((ll)floor(P1[i].y+0.49)-a1b2)%mod;
    a2b2=((ll)floor(P2[i].x+0.49)-a1b1)%mod;
    f[i]=((a1b1*lim+(a1b2+a2b1))*lim+a2b2)%mod;
    f[i]=(f[i]+mod)%mod;
  }for (int i=0;i<len;i++)
    P1[i]=P2[i]=Q[i]=(CP){0.0,0.0};
}
long long C(int n,int m)
{return fac[n]*inv[m]%mod*inv[n-m]%mod;}
void Init(int limit)
{
  inv[1]=inv[0]=fac[0]=1;
  for (int i=1;i<=limit;i++)fac[i]=fac[i-1]*i%mod;
  for (int i=2;i<=limit;i++)
    inv[i]=inv[mod%i]*(mod-mod/i)%mod;
  for (int i=2;i<=limit;i++)inv[i]=inv[i-1]*inv[i]%mod;
}
ll p[Maxn<<2];
void rev(ll *f,int len)
{
  for (int i=0;i<len;i++)p[i]=f[i];
  for (int i=0;i<len;i++)f[len-i-1]=p[i];
}
//求出F(x-c) 
void fminus(ll *s,ll *f,int len,int c)
{
  c=mod-c;
  for (int i=0;i<len;i++)
    p[len-i-1]=f[i]*fac[i]%mod;
  ll buf=1;
  for (int i=0;i<len;i++,buf=buf*c%mod)
    s[i]=buf*inv[i]%mod;
  times(p,s,len,len);
  for (int i=0;i<len;i++)s[len-i-1]=p[i]*inv[len-i-1]%mod;
  for (int i=len;i<len+len;i++)s[i]=0;
}
ll f[Maxn<<2],s[Maxn<<2];
void solve(ll *f,int n)
{
  if (n==1){f[0]=0;f[1]=1;}
  else if (n&1){
    solve(f,n-1);f[n]=0;
    //再乘上(x-n+1)就好了
    for (int i=n;i>0;i--)
      f[i]=(f[i-1]+(mod-n+1)*f[i])%mod;
    f[0]=f[0]*(mod-n+1)%mod;
  }else {
    solve(f,n/2);
    //S(x)=F(x+n/2)
    fminus(s,f,n/2+1,n/2);
    times(f,s,n/2+1,n+1);
  }
}
void invp(ll *f,int len)
{
  for (int i=0;i<k+1;i++)s[i]=p[i]=0;
  //注意清空 
  ll *r=s,*rr=p;
  int n=1;for(;n<len;n<<=1);
  rr[0]=powM(f[0]);
  for (int len=2;len<=n;len<<=1){
    for (int i=0;i<len;i++)
      r[i]=rr[i]*2%mod;
    times(rr,rr,len/2,len);
    times(rr,f,len,len);
    for (int i=0;i<len;i++)
      rr[i]=(r[i]-rr[i]+mod)%mod;
  }for (int i=0;i<len;i++)
    f[i]=rr[i];
}
long long ans,sum;
int main()
{
  scanf("%d%d",&n,&k);
  if (k>n){printf("0");return 0;}
  for (int i=1,sav;i<=n;i++){
    scanf("%d",&sav);
    ans=(ans+sav)%mod;
  }k--;
  Init(n);solve(f,k+1);
  for (int i=0;i<k+1;i++)f[i]=f[i+1];
  rev(f,k+1);
  for (int i=n-k+1;i<k+1;i++)f[i]=0;
  for (int i=k+1;i<n-k+1;i++)f[i]=0;
  invp(f,n-k+1);
  for (int i=n+1;i>=k;i--)f[i]=f[i-k];
  for (int i=0;i<k;i++)f[i]=0;
  for (int i=1;i<n+1;i++)
    sum=(sum+1ll*i*C(n-1,i-1)%mod*f[n-i])%mod;
  printf("%lld",sum*ans%mod);
  return 0;
}

另一种方法是考虑大力推式,直到我们逃离卷积。

p=\sum\limits_{i=1}^ni\dbinom{n-1}{i-1}\begin{Bmatrix}n-i\\k-1\end{Bmatrix}

回顾 : m^n=\sum\limits_{j=0}^m\dbinom{m}{j}\begin{Bmatrix}n\\j\end{Bmatrix}j!
\Longrightarrow\begin{Bmatrix}n\\m\end{Bmatrix}=\frac{1}{m!}\sum\limits_{j=0}^m\dbinom{m}{j}(-1)^{m-j}i^n

=\sum\limits_{i=1}^ni\dbinom{n-1}{i-1}\frac{1}{(k-1)!}\sum\limits_{j=0}^{k-1}\dbinom{k-1}{j}(-1)^{k-1-j}j^{n-i}

拆开组合数,约掉一些东西。

=\sum\limits_{i=1}^ni\dbinom{n-1}{i-1}\sum\limits_{j=0}^{k-1}\dfrac{(-1)^{k-1-j}j^{n-i}}{j!(k-1-j)!}

交换和式可得 :

=\sum\limits_{j=0}^{k-1}\dfrac{(-1)^{k-1+j}}{j!(k-1-j)!}\sum\limits_{i=1}^ni\dbinom{n-1}{i-1}j^{n-i}

现在来单独关注\sum\limits_{i=1}^ni\binom{n-1}{i-1}j^{n-i}
- 回顾 : \dbinom{n}{k}=\dfrac{n}{k}\dbinom{n-1}{k-1}
注意到有个常数i阻碍了处理,而又不能直接吸收掉,先拆成1+(i-1)
=\sum\limits_{i=1}^n\dbinom{n-1}{i-1}j^{n-i}+\sum\limits_{i=1}^n(i-1)\dbinom{n-1}{i-1}j^{n-i}
在后半部分作吸收,得到
=\sum\limits_{i=1}^n\dbinom{n-1}{i-1}j^{n-i}+(n-1)\sum\limits_{i=1}^n\dbinom{n-2}{i-2}j^{n-i} =\sum\limits_{i=0}^{n-1}\dbinom{n-1}{i}j^{n-i-1}+(n-1)\sum\limits_{i=0}^{n-2}\dbinom{n-2}{i}j^{n-i-2} =\sum\limits_{i=0}^{n-1}\dbinom{n-1}{i}j^{(n-1)-i}+(n-1)\sum\limits_{i=0}^{n-2}\dbinom{n-2}{i}j^{(n-2)-i}
现在可以使用二项式定理。
=(j+1)^{n-1}+(n-1)(j+1)^{n-2}

然后随便怎么做就O(n\log n)的了。

还有一种更妙的方法 : 愿组合意义指引你!

我们观察 W(S)=|S|\sum\limits_{x\in S}w_x 前面的系数|S|有点诡异,尝试用组合意义拼凑。

事实上,可以视作划分出的集合中,每两对点之间会贡献,贡献是两者价值和。

现在,我们对于某个元素,考虑所有元素与其配对的贡献系数。

自己和自己 : \begin{Bmatrix}n\\k\end{Bmatrix} (所有方案中都能与自己配对)

其他点 : (n-1)\begin{Bmatrix}n-1\\k\end{Bmatrix} (先分好集合,再放进去)

于是p=\begin{Bmatrix}n\\k\end{Bmatrix}+(n-1)\begin{Bmatrix}n-1\\k\end{Bmatrix}

```cpp ```