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 5Sample Output
67584000
开局一个式子, 装备全靠混凝土
首先我们发现贡献可以非常容易地根据每个点的情况分开来算, 我们有:
$$
F(n)=n\sum_{i=1}^{n-1}{n-1 \choose i}i^k2^{{n\choose 2}-(n-1)}
$$
含义就是枚举某个点的度数, 算出和它相关的连边方案, 其他边随便连. 点和点之间等价, 直接乘 $n$.
接着我们发现和式内部有些和 $i$ 无关的东西, 我们把它扥出来:
$$
F(n)=n2^{{n\choose 2}-(n-1)}\sum_{i=1}^{n-1}{n-1\choose i}i^k
$$
前面部分显然非常好求, 我们把重点放在后面那个和式上, 设:
$$
G(n)=\sum_{i=1}^n{n\choose i}i^k
$$
这个求和式的项数是和 $n$ 有关的, 然而 $n$ 巨大无比根本跑不出来...只有一个 $k$ 的数据范围还有救...
我们考虑如何把求和上界向 $k$ 的方向转化.
翻阅混凝土数学, 我们找到了第二类Stirling反演公式 $(6.10)$:
$$
x^k=\sum_{i=0}^k\begin{Bmatrix}k\\i\end{Bmatrix}x^{\underline i}
$$
扔进去试试看:
$$
G(n)=\sum_{i=1}^n{n\choose i}\sum_{r=0}^k\begin{Bmatrix}k\\r\end{Bmatrix}i^{\underline r}
$$
按照求和式套路把所有东西放在最里层, 对外层求和指标乱搞再把和里层指标无关的东西抻回来:
$$
\begin{aligned}
G(n)&=\sum_{i=1}^n\sum_{r=0}^k{n\choose i}\begin{Bmatrix}k\\r\end{Bmatrix}i^{\underline r} \\
&=\sum_{r=0}^k\sum_{i=1}^n{n\choose i}\begin{Bmatrix}k\\r\end{Bmatrix}i^{\underline r} \\
&=\sum_{r=0}^k\begin{Bmatrix}k\\r\end{Bmatrix}\sum_{i=1}^n{n\choose i}i^{\underline r}
\end{aligned}
$$
我们看那个下降阶乘幂在一个二项式系数旁边就非常不爽, 我们把它换成二项式系数来方便比对公式:
$$
G(n)=\sum_{r=0}^k\begin{Bmatrix}k\\r\end{Bmatrix}\sum_{i=1}^n{n\choose i}{i\choose r}r!
$$
然后我们果然在混凝土中匹配到了一个公式 $(5.21)$!
$$
{r\choose m}{m\choose k}={r\choose k}{r-k\choose m-k}
$$
代进去化一化:
$$
\begin{aligned}
G(n)&=\sum_{r=0}^k\begin{Bmatrix}k\\r\end{Bmatrix}\sum_{i=1}^n{n\choose r}{n-r\choose i-r}r!\\
&=\sum_{r=0}^k\begin{Bmatrix}k\\r\end{Bmatrix}r!{n\choose r}\sum_{i=1}^n{n-r\choose i-r}\\
&=\sum_{r=0}^k\begin{Bmatrix}k\\r\end{Bmatrix}r!{n\choose r}2^{n-r}
\end{aligned}
$$
完美!
现在问题是怎么求 $\begin{Bmatrix}k\\r\end{Bmatrix}$. 它可能并没有什么封闭形式.
回想 $\begin{Bmatrix}k\\r\end{Bmatrix}$ 到底表示什么:
符号 $\begin{Bmatrix}n\\k\end{Bmatrix}$ 表示将一个有 $n$ 个物品的集合划分成 $k$ 个非空子集的方案数
如果没有集合非空的限制而且集合带有 1~k 的标号, 那么我们可以发现答案就是 $k^n$(给 $n$ 个元素分配一个集合标号)
枚举有多少集合是空的, 容斥一下并消序得到:
$$
\begin{Bmatrix}n\\k\end{Bmatrix}=\frac 1 {k!}\sum_{i=0}^k{k\choose i}(-1)^i(k-i)^n
$$
我们发现后面的和式就是数列 $\left \langle (-1)^k \right \rangle$ 和 $\left \langle k^n \right \rangle$ 的二项卷积, 法法塔一发即可 $O(n\log n)$ 求出所有 $\begin{Bmatrix} n \\ k \end{Bmatrix}$.
求完Stirling数回代即可. 复杂度 $O(k\log 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\n",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 }
原文:https://www.cnblogs.com/rvalue/p/10240700.html