题解

思路概述

• 题目要求计算最小编号黑球 A 的期望，并给出多项式 F(A) 的取值。设两次随机着色的分布转化为“最左侧黑球在位置 x”的概率 prob[x]，则答案是 Σ F(x)·prob[x]。
• 为了高效求出这组概率，代码使用组合数与容斥公式，将“前 x-1 个位置均未被涂黑、位置 x 被涂黑”转化为若干组合项，并借助卷积与 NTT 进行快速计算。
• 多项式系数都是模 998244353 的值；程序先对输入 (F(0)…F(m)) 进行差分求出多项式系数，再通过两次 NTT 做卷积，把公式化为多项式乘法；最后根据公式累加求出 Σ F(x)·prob[x] 并乘以题目指定的组合系数。

复杂度

• 归约到一次 O(m log m) 的 NTT 卷积及若干线性遍历，足以应对 m ≤ 10^6。
• 代码中 init()、DIT/DIF 等函数实现了模数 998244353 下的 NTT 变换。

#include<map>
#include<set>
#include<ctime>
#include<cmath>
#include<queue>
#include<bitset>
#include<cstdio>
#include<vector>
#include<random>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#define ll long long
using namespace std;
#define I ll
#define her1 20081214
#define IV void
#define cht 998244353
#define ld double
#define ull unsigned long long
#define cp(x,y)memcpy(x,y,sizeof y)
#define mem(x,val)memset(x,val,sizeof x)
#define D(i,j,n)for(register int i=j;i>=n;i--)
#define E(i,now)for(register int i=first[now],v=e[i].to;i;v=e[i=e[i].nxt].to)
#define F(i,j,n)for(register int i=j;i<=n;i++)
#define DL(i,j,n)for(register i64 i=j;i>=n;i--)
#define EL(i,now)for(register i64 i=first[now];i;i=e[i].nxt)
#define FL(i,j,n)for(register i64 i=j;i<=n;i++)
//#define D(i,j,n)for(int i=j;i>=n;i--)
//#define E(i,now)for(int i=first[now];i;i=e[i].nxt)
//#define F(i,j,n)for(int i=j;i<=n;i++)
//#define DL(i,j,n)for(register ll i=j;i>=n;i--)
//#define EL(i,now)for(register ll i=first[now];i;i=e[i].nxt)
//#define FL(i,j,n)for(register ll i=j;i<=n;i++)
ll read(){
	ll ans=0,f=1;
	char c=getchar();
	while(c<'0'||c>'9'){
		if(c=='-')f=-1;
		c=getchar();
	}
	while(c>='0'&&c<='9')ans=ans*10+c-'0',c=getchar();
	return ans*f;
}
#undef ll
#include "assert.h"
mt19937_64 rnd(her1);
#include "functional"
using i64 = long long;

const int MAX = (1<<22);
const int maxn = (1<<22)+5;
const int maxm = (1<<22)+5;
#define in(i,V) F(i,0,(i64)V.size()-1)
#define My_assert(expr,tips) ((expr)?(void(0)):(puts(tips),exit(0)))
IV cadd(i64&x,i64 val){x=(x+val)%cht;}

i64 w[maxm],fac[maxn],ifac[maxn];
i64 qpow(i64 n,i64 base=cht-2){
	i64 ans=1;
	while(base){
		if(base&1)ans=ans*n%cht;
		n=n*n%cht;base>>=1;
	}
	return ans;
}
IV init(i64 limit){
	for(int i=1,j,k;i<limit;i<<=1)
		for(w[j=i]=1,k=qpow(3,(cht-1)/(i<<1)),j++;j<(i<<1);j++)
			w[j]=w[j-1]*k%cht;
}
IV DIT(i64*a,i64 limit){
	for(i64 i,j,k=limit>>1,L,*W,*x,*y,z;k;k>>=1)
		for(L=k<<1,i=0;i<limit;i+=L)for(j=0,W=w+k,x=a+i,y=x+k;j<k;j++,W++,x++,y++)
			*y=(*x+cht-(z=*y))* *W%cht,*x=(*x+z)%cht;
}
IV DIF(i64*a,i64 limit){
	for(i64 i,j,k=1,L,*W,*x,*y,z;k<limit;k<<=1)
		for(L=k<<1,i=0;i<limit;i+=L)for(j=0,W=w+k,x=a+i,y=x+k;j<k;j++,W++,x++,y++)
			z=1ll* *W* *y%cht,*y=(*x-z)%cht,*x=(*x+z)%cht;
	reverse(a+1,a+limit);
	i64 inv=qpow(limit);
	F(i,0,limit-1)a[i]=a[i]*inv%cht;
}
#define poly vector<i64>

poly operator+(const poly&A,const poly&B){
	poly C;C.resize(max(A.size(),B.size()));
	in(i,C)C[i]=0;
	in(i,A)cadd(C[i],A[i]);
	in(i,B)cadd(C[i],B[i]);
	return C;
}
poly operator-(const poly&A,const poly&B){
	poly C;C.resize(max(A.size(),B.size()));
	in(i,C)C[i]=0;
	in(i,A)cadd(C[i],A[i]);
	in(i,B)cadd(C[i],-B[i]);
	return C;
}
poly operator*(const poly&A,const poly&B){
	static i64 f[maxm],g[maxm];
	if(A.empty()||B.empty())return{};
	My_assert(A.size()&&B.size(),"multiply with a empty poly.");
	i64 limit=1;
	while(limit<=A.size()+B.size())
		limit<<=1;
	F(i,0,limit-1)f[i]=g[i]=0;
	in(i,A)f[i]=A[i];in(i,B)g[i]=B[i];
	DIT(f,limit);DIT(g,limit);
	F(i,0,limit-1)f[i]=f[i]*g[i]%cht;DIF(f,limit);
	poly C;C.resize(A.size()+B.size()-1);
	in(i,C)C[i]=f[i];
	return C;
}
poly operator*(const poly&A,const i64&k){
	poly Ans=A;for(auto&x:Ans)x=x*k%cht;
	return Ans;
}
#define vi vector<i64>
IV print(poly A){for(auto x:A)cout<<x<<' ';puts("");}
IV remod(poly&A,i64 n){while(A.size()>n)A.pop_back();}
IV init(){
	fac[0]=1;F(i,1,MAX)fac[i]=fac[i-1]*i%cht;
	ifac[MAX]=qpow(fac[MAX]);D(i,MAX-1,0)ifac[i]=ifac[i+1]*(i+1)%cht;
}
i64 C(i64 n,i64 m){
	if(n<m||n<0||m<0)return 0;
	return fac[n]*ifac[m]%cht*ifac[n-m]%cht;
}
i64 n,m,f[maxn],g[maxn];
i64 calc(i64 x){
	i64 mul=x,sum=0;
	if(x<=m)return f[x];
	F(i,0,m){
		i64 M=f[i];
		F(j,0,m)if(i!=j)M=M*(x-j)%cht*qpow(i-j)%cht;
		sum=(sum+M)%cht;
	}
	return sum;
}
i64 sqr(i64 x){return x*x%cht;}

i64 B[maxn];
poly IDFT(poly A,i64 n){
	poly B(ifac,ifac+n);
	in(i,A)A[i]=A[i]*ifac[i]%cht;
	F(i,0,n-1)if(i&1)B[i]=-B[i];
	A=A*B;remod(A,n);return A;
}
poly DFT(poly A,i64 n){
	poly B(ifac,ifac+n);A=A*B;remod(A,n);
	F(i,0,n-1)A[i]=A[i]*fac[i]%cht;return A;
}
i64 pre[maxn],inv[maxn],A[maxn];
int main(){
	// freopen("1.in","r",stdin);
	// freopen("1.out","w",stdout);
	init(1<<21);init();

	n=read();m=read();
	poly A(m+1),B(ifac,ifac+m+1);
	F(i,0,m)A[i]=read()*ifac[i]%cht;
	F(i,0,m)if(i&1)B[i]=-B[i];A=A*B;
	F(i,0,m)A[i]=A[i]*fac[i]%cht;remod(A,m+1);
	poly T(m+1);
	F(k,0,m)T[k]=C(k+m,m)*C(m,m-k)%cht;
	A=A*T;i64 Ans=0,mul=1;
	F(s,0,3*m){
		if(s>=m)cadd(Ans,mul*ifac[s]%cht*A[s-m]);
		mul=mul*(n-s)%cht;
	}
	cout<<(Ans+cht)%cht;
	return 0;
}