zoukankan      html  css  js  c++  java
  • [BZOJ 5093] 图的价值

    5093: [Lydsy1711月赛]图的价值

    Time Limit: 30 Sec  Memory Limit: 256 MB
    Submit: 544  Solved: 273
    [Submit][Status][Discuss]

    Description

    “简单无向图”是指无重边、无自环的无向图(不一定连通)。
    一个带标号的图的价值定义为每个点度数的k次方的和。
    给定n和k,请计算所有n个点的带标号的简单无向图的价值之和。
    因为答案很大,请对998244353取模输出。
     

    Input

    第一行包含两个正整数n,k(1<=n<=10^9,1<=k<=200000)。
     

    Output

     输出一行一个整数,即答案对998244353取模的结果。

     

    Sample Input

    6 5

    Sample Output

    67584000

    开局一个式子, 装备全靠混凝土

    首先我们发现贡献可以非常容易地根据每个点的情况分开来算, 我们有:
    $$
    F(n)=nsum_{i=1}^{n-1}{n-1 choose i}i^k2^{{nchoose 2}-(n-1)}
    $$
    含义就是枚举某个点的度数, 算出和它相关的连边方案, 其他边随便连. 点和点之间等价, 直接乘 $n$.

    接着我们发现和式内部有些和 $i$ 无关的东西, 我们把它扥出来:
    $$
    F(n)=n2^{{nchoose 2}-(n-1)}sum_{i=1}^{n-1}{n-1choose i}i^k
    $$
    前面部分显然非常好求, 我们把重点放在后面那个和式上, 设:
    $$
    G(n)=sum_{i=1}^n{nchoose i}i^k
    $$
    这个求和式的项数是和 $n$ 有关的, 然而 $n$ 巨大无比根本跑不出来...只有一个 $k$ 的数据范围还有救...

    我们考虑如何把求和上界向 $k$ 的方向转化.

    翻阅混凝土数学, 我们找到了第二类Stirling反演公式 $(6.10)$:
    $$
    x^k=sum_{i=0}^kegin{Bmatrix}k\iend{Bmatrix}x^{underline i}
    $$
    扔进去试试看:
    $$
    G(n)=sum_{i=1}^n{nchoose i}sum_{r=0}^kegin{Bmatrix}k\rend{Bmatrix}i^{underline r}
    $$
    按照求和式套路把所有东西放在最里层, 对外层求和指标乱搞再把和里层指标无关的东西抻回来:
    $$
    egin{aligned}
    G(n)&=sum_{i=1}^nsum_{r=0}^k{nchoose i}egin{Bmatrix}k\rend{Bmatrix}i^{underline r} \
    &=sum_{r=0}^ksum_{i=1}^n{nchoose i}egin{Bmatrix}k\rend{Bmatrix}i^{underline r} \
    &=sum_{r=0}^kegin{Bmatrix}k\rend{Bmatrix}sum_{i=1}^n{nchoose i}i^{underline r}
    end{aligned}
    $$
    我们看那个下降阶乘幂在一个二项式系数旁边就非常不爽, 我们把它换成二项式系数来方便比对公式:
    $$
    G(n)=sum_{r=0}^kegin{Bmatrix}k\rend{Bmatrix}sum_{i=1}^n{nchoose i}{ichoose r}r!
    $$
    然后我们果然在混凝土中匹配到了一个公式 $(5.21)$!
    $$
    {rchoose m}{mchoose k}={rchoose k}{r-kchoose m-k}
    $$
    代进去化一化:
    $$
    egin{aligned}
    G(n)&=sum_{r=0}^kegin{Bmatrix}k\rend{Bmatrix}sum_{i=1}^n{nchoose r}{n-rchoose i-r}r!\
    &=sum_{r=0}^kegin{Bmatrix}k\rend{Bmatrix}r!{nchoose r}sum_{i=1}^n{n-rchoose i-r}\
    &=sum_{r=0}^kegin{Bmatrix}k\rend{Bmatrix}r!{nchoose r}2^{n-r}
    end{aligned}
    $$
    完美!

    现在问题是怎么求 $egin{Bmatrix}k\rend{Bmatrix}$. 它可能并没有什么封闭形式.

    回想 $egin{Bmatrix}k\rend{Bmatrix}$ 到底表示什么:

    符号 $egin{Bmatrix}n\kend{Bmatrix}$ 表示将一个有 $n$ 个物品的集合划分成 $k$ 个非空子集的方案数

    如果没有集合非空的限制而且集合带有 1~k 的标号, 那么我们可以发现答案就是 $k^n$(给 $n$ 个元素分配一个集合标号)

    枚举有多少集合是空的, 容斥一下并消序得到:
    $$
    egin{Bmatrix}n\kend{Bmatrix}=frac 1 {k!}sum_{i=0}^k{kchoose i}(-1)^i(k-i)^n
    $$
    我们发现后面的和式就是数列 $left langle (-1)^k ight angle$ 和 $left langle k^n ight angle$ 的二项卷积, 法法塔一发即可 $O(nlog n)$ 求出所有 $egin{Bmatrix} n \ k end{Bmatrix}$.

    求完Stirling数回代即可. 复杂度 $O(klog k)$.

    代码实现

    就用了个多项式乘法还套全家桶真是吃枣药丸

      1 #include <bits/stdc++.h>
      2 
      3 const int G=3;
      4 const int DFT=1;
      5 const int IDFT=-1;
      6 const int MAXN=550000;
      7 const int MOD=998244353;
      8 const int INV2=(MOD+1)>>1;
      9 const int PHI=MOD-1;
     10 
     11 typedef std::vector<int> Poly;
     12 
     13 Poly Sqrt(Poly);
     14 void Read(Poly&);
     15 Poly Inverse(Poly);
     16 Poly Ln(const Poly&);
     17 Poly Exp(const Poly&);
     18 void Print(const Poly&);
     19 void NTT(Poly&,int,int);
     20 Poly Pow(const Poly&,int);
     21 Poly Integral(const Poly&);
     22 Poly Derivative(const Poly&);
     23 Poly operator*(Poly,Poly);
     24 Poly operator/(Poly,Poly);
     25 Poly operator%(Poly,Poly);
     26 Poly operator+(const Poly&,const Poly&);
     27 Poly operator-(const Poly&,const Poly&);
     28 
     29 int Cn[MAXN];
     30 int rev[MAXN];
     31 int inv[MAXN];
     32 int fact[MAXN];
     33 
     34 int NTTPre(int);
     35 int Sqrt(int,int);
     36 int Pow(int,int,int);
     37 int Log(int,int,int);
     38 int ExGCD(int,int,int&,int&);
     39 
     40 int main(){
     41     int n,k;
     42     scanf("%d%d",&n,&k);
     43     int ans=1ll*n*Pow(2,1ll*(n-1)*(n-2)/2%(MOD-1),MOD)%MOD;
     44     --n;
     45     Poly a(k+1),b(k+1);
     46     fact[0]=a[0]=Cn[0]=1;
     47     for(int i=1;i<=k;i++){
     48         fact[i]=1ll*fact[i-1]*i%MOD;
     49         inv[i]=Pow(fact[i],MOD-2,MOD);
     50         a[i]=1ll*(i&1?MOD-1:1)*inv[i]%MOD;
     51         b[i]=1ll*Pow(i,k,MOD)*inv[i]%MOD;
     52         Cn[i]=1ll*Cn[i-1]*(n-i+1)%MOD*Pow(i,MOD-2,MOD)%MOD;
     53     }
     54     Poly s=a*b;
     55     int g=0;
     56     for(int i=0;i<=k;i++)
     57         (g+=1ll*s[i]*fact[i]%MOD*Cn[i]%MOD*Pow(2,n-i,MOD)%MOD)%=MOD;
     58     ans=1ll*g*ans%MOD;
     59     printf("%d
    ",ans);
     60     return 0;
     61 }
     62 
     63 void Read(Poly& a){
     64     for(auto& i:a)
     65         scanf("%d",&i);
     66 }
     67 
     68 void Print(const Poly& a){
     69     for(auto i:a)
     70         printf("%d ",i);
     71     puts("");
     72 }
     73 
     74 Poly Pow(const Poly& a,int k){
     75     Poly log=Ln(a);
     76     for(auto& i:log)
     77         i=1ll*i*k%MOD;
     78     return Exp(log);
     79 }
     80 
     81 Poly Sqrt(Poly a){
     82     int len=a.size();
     83     if(len==1)
     84         return Poly(1,Sqrt(a[0],MOD));
     85     else{
     86         Poly b=a;
     87         b.resize((len+1)>>1);
     88         b=Sqrt(b);
     89         b.resize(len);
     90         Poly inv=Inverse(b);
     91         int bln=NTTPre(inv.size()+a.size());
     92         NTT(a,bln,DFT);
     93         NTT(inv,bln,DFT);
     94         for(int i=0;i<bln;i++)
     95             a[i]=1ll*a[i]*INV2%MOD*inv[i]%MOD;
     96         NTT(a,bln,IDFT);
     97         for(int i=0;i<len;i++)
     98             b[i]=(1ll*b[i]*INV2%MOD+a[i])%MOD;
     99         return b;
    100     }
    101 }
    102 
    103 Poly Exp(const Poly& a){
    104     size_t len=1;
    105     Poly ans(1,1),one(1,1);
    106     while(len<(a.size()<<1)){
    107         len<<=1;
    108         Poly b=a;
    109         b.resize(len);
    110         ans=ans*(one-Ln(ans)+b);
    111         ans.resize(len);
    112     }
    113     ans.resize(a.size());
    114     return ans;
    115 }
    116 
    117 Poly Ln(const Poly& a){
    118     Poly ans=Integral(Derivative(a)*Inverse(a));
    119     ans.resize(a.size());
    120     return ans;
    121 }
    122 
    123 Poly Integral(const Poly& a){
    124     int len=a.size();
    125     Poly ans(len+1);
    126     for(int i=1;i<len;i++)
    127         ans[i]=1ll*a[i-1]*Pow(i,MOD-2,MOD)%MOD;
    128     return ans;
    129 }
    130 
    131 Poly Derivative(const Poly& a){
    132     int len=a.size();
    133     Poly ans(len-1);
    134     for(int i=1;i<len;i++)
    135         ans[i-1]=1ll*a[i]*i%MOD;
    136     return ans;
    137 }
    138 
    139 Poly operator/(Poly a,Poly b){
    140     int n=a.size()-1,m=b.size()-1;
    141     Poly ans(1);
    142     if(n>=m){
    143         std::reverse(a.begin(),a.end());
    144         std::reverse(b.begin(),b.end());
    145         b.resize(n-m+1);
    146         ans=Inverse(b)*a;
    147         ans.resize(n-m+1);
    148         std::reverse(ans.begin(),ans.end());
    149     }
    150     return ans;
    151 }
    152 
    153 Poly operator%(Poly a,Poly b){
    154     int n=a.size()-1,m=b.size()-1;
    155     Poly ans;
    156     if(n<m)
    157         ans=a;
    158     else
    159         ans=a-(a/b)*b;
    160     ans.resize(m);
    161     return ans;
    162 }
    163 
    164 Poly operator*(Poly a,Poly b){
    165     int len=a.size()+b.size()-1;
    166     int bln=NTTPre(len);
    167     NTT(a,bln,DFT);
    168     NTT(b,bln,DFT);
    169     for(int i=0;i<bln;i++)
    170         a[i]=1ll*a[i]*b[i]%MOD;
    171     NTT(a,bln,IDFT);
    172     a.resize(len);
    173     return a;
    174 }
    175 
    176 Poly operator+(const Poly& a,const Poly& b){
    177     Poly ans(std::max(a.size(),b.size()));
    178     std::copy(a.begin(),a.end(),ans.begin());
    179     for(size_t i=0;i<b.size();i++)
    180         ans[i]=(ans[i]+b[i])%MOD;
    181     return ans;
    182 }
    183 
    184 Poly operator-(const Poly& a,const Poly& b){
    185     Poly ans(std::max(a.size(),b.size()));
    186     std::copy(a.begin(),a.end(),ans.begin());
    187     for(size_t i=0;i<b.size();i++)
    188         ans[i]=(ans[i]+MOD-b[i])%MOD;
    189     return ans;
    190 }
    191 
    192 Poly Inverse(Poly a){
    193     int len=a.size();
    194     if(len==1)
    195         return Poly(1,Pow(a[0],MOD-2,MOD));
    196     else{
    197         Poly b(a);
    198         b.resize((len+1)>>1);
    199         b=Inverse(b);
    200         int bln=NTTPre(b.size()*2+a.size());
    201         NTT(a,bln,DFT);
    202         NTT(b,bln,DFT);
    203         for(int i=0;i<bln;i++)
    204             b[i]=(2ll*b[i]%MOD-1ll*b[i]*b[i]%MOD*a[i]%MOD+MOD)%MOD;
    205         NTT(b,bln,IDFT);
    206         b.resize(len);
    207         return b;
    208     }
    209 }
    210 
    211 void NTT(Poly& a,int len,int opt){
    212     a.resize(len);
    213     for(int i=0;i<len;i++)
    214         if(rev[i]>i)
    215             std::swap(a[i],a[rev[i]]);
    216     for(int i=1;i<len;i<<=1){
    217         int step=i<<1;
    218         int wn=Pow(G,(PHI+opt*PHI/step)%PHI,MOD);
    219         for(int j=0;j<len;j+=step){
    220             int w=1;
    221             for(int k=0;k<i;k++,w=1ll*w*wn%MOD){
    222                 int x=a[j+k];
    223                 int y=1ll*a[j+k+i]*w%MOD;
    224                 a[j+k]=(x+y)%MOD;
    225                 a[j+k+i]=(x-y+MOD)%MOD;
    226             }
    227         }
    228     }
    229     if(opt==IDFT){
    230         int inv=Pow(len,MOD-2,MOD);
    231         for(int i=0;i<len;i++)
    232             a[i]=1ll*a[i]*inv%MOD;
    233     }
    234 }
    235 
    236 int NTTPre(int n){
    237     int bln=1,bct=0;
    238     while(bln<n){
    239         bln<<=1;
    240         ++bct;
    241     }
    242     for(int i=0;i<bln;i++)
    243         rev[i]=(rev[i>>1]>>1)|((i&1)<<(bct-1));
    244     return bln;
    245 }
    246 
    247 inline int Pow(int a,int n,int p){
    248     int ans=1;
    249     while(n>0){
    250         if(n&1)
    251             ans=1ll*a*ans%p;
    252         a=1ll*a*a%p;
    253         n>>=1;
    254     }
    255     return ans;
    256 }
    257 
    258 int ExGCD(int a,int b,int& x,int& y){
    259     if(b==0){
    260         x=1,y=0;
    261         return a;
    262     }
    263     else{
    264         int gcd=ExGCD(b,a%b,y,x);
    265         y-=x*(a/b);
    266         return gcd;
    267     }
    268 }
    269 
    270 inline int Sqrt(int a,int p){
    271     int s=Pow(G,Log(G,a,p)>>1,p);
    272     return std::min(s,MOD-s);
    273 }
    274 
    275 inline int Log(int a,int x,int p){
    276     int s=sqrt(p+0.5);
    277     int inv=Pow(Pow(a,s,p),p-2,p);
    278     std::unordered_map<int,int> m;
    279     m[1]=0;
    280     int pow=1;
    281     for(int i=1;i<s;i++){
    282         pow=1ll*a*pow%p;
    283         if(!m.count(pow))
    284             m[pow]=i;
    285     }
    286     for(int i=0;i<s;i++){
    287         if(m.count(x))
    288             return i*s+m[x];
    289         x=1ll*x*inv%MOD;
    290     }
    291     return -1;
    292 }
    BZOJ 5093

    日常图包

  • 相关阅读:
    读书笔记——吴军《态度》
    JZYZOJ1237 教授的测试 dfs
    NOI1999 JZYZOJ1289 棋盘分割 dp 方差的数学结论
    [JZYZOJ 1288][洛谷 1005] NOIP2007 矩阵取数 dp 高精度
    POJ 3904 JZYZOJ 1202 Sky Code 莫比乌斯反演 组合数
    POJ2157 Check the difficulty of problems 概率DP
    HDU3853 LOOPS 期望DP 简单
    Codeforces 148D. Bag of mice 概率dp
    POJ3071 Football 概率DP 简单
    HDU4405 Aeroplane chess 飞行棋 期望dp 简单
  • 原文地址:https://www.cnblogs.com/rvalue/p/10240700.html
Copyright © 2011-2022 走看看