【模拟赛12.22】Graph

Time Limit 3s
Memory Limit 512M
Description
认定两个黑白二分图不同当且仅当某两个点在一个二分图中有边在另一个中没边，或者在两个图中边颜色不同，求n个点的不同的黑白二分图的个数
Input
第一行T数据组数，然后T行每行n
Output
每行一个答案对998244353取模
Sample

Input
3
2
3
10
Output
3
19
435363530

Data Constraint
$n\le 10^5$

---

考虑将n分为多个联通块，设一个大小为n的联通快的方案数的生成函数为$f(x)$,$EGF(f)=F(x)$
答案就是$n!\ast [x^n]e^{F(x)}$
这里$n!$乘的是里面的$EGF$，而外面$e^{F(x)}=\sum \frac{F(x{)}^i}{i!}$，除$i!$是为了避免选取的$i$个被重复计算$i!$次
显然，不好推
考虑$g(x)$为${\sum }_{x+y=n}\left(\genfrac{}{}{0ex}{}{n}{x}\right)3^{xy}$的生成函数，其EGF为$G(x)$，这样是不考虑这一块是否联通的，那么如果这里有k个联通快，由于左右两部可以互换，重复算了$2^k$次，有关系式$G(x)=\sum \frac{F(x{)}^i{2}^i}{i!}=e^{2F(x)}$，那么$e^{F(x)}=\sqrt{G(x)}$

---

Part 1
如何求G(x)呢
二次剩余？？很遗憾3在摸998244353下无二次剩余。。。
$$\begin{array}{rl}[x^n]G(x)& =\frac{{\sum }_{i+j=n}\left(\genfrac{}{}{0ex}{}{n}{i}\right)3^{ij}}{n!}\\ & =3^{\frac{n^2}{2}}\ast \sum _{i+j=n}\frac{3^{\frac{-i^2}{2}}}{i!}\times \frac{3^{\frac{-j^2}{2}}}{j!}\end{array}$$这个分奇偶数讨论一下就好了

---

Part 2
多项式开根
假设原式为A(x)
假设
$G(x{)}^2\equiv A(x)(\mathrm{m}\mathrm{o}\mathrm{d}{x}^{\frac{n}{2}})$
$B(x{)}^2\equiv A(x)(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)$
这时
$$\begin{array}{rl}B(x)-G(x)& \equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^{\frac{n}{2}})\\ (B(x)-G(x){)}^2& \equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)\\ B^2(x)-2B(x)G(x)+G^2(x)& \equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)\\ A(x)-2B(x)G(x)+G^2(x)& \equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)\\ B(x)& \equiv \frac{G(x)}{2}+\frac{A(x)}{2G(x)}(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)\end{array}$$
这时就只需要知道答案的常数项就好了
这题A的常数项为1，比较特殊
不然用二次剩余具体是这样：
勒让德符号$(\frac{n}{p})$
$$(\frac{n}{p})=n^{\frac{p-1}{2}}=\{\begin{array}{ll}1& ,存在二次剩余\\ 0& ,n=0\\ -1& ,不存在二次剩余\end{array}(\mathrm{m}\mathrm{o}\mathrm{d}p)$$
加入n有实根
要找$a^2=n(\mathrm{m}\mathrm{o}\mathrm{d}p)$，先随机一个$(\frac{t^2-n}{p})=-1$，令$w=\sqrt{t^2-n}(可以看作是复数根)$，那么$a=(1-w{)}^{\frac{p+1}{2}}(显然这个只有实部)$

---

Part 3
如何去求$\frac{A(x)}{2G(x)}$？考虑求$\frac{1}{G(x)}$
多项式求逆！
同开根的思路
设
$B(x)A(x)\equiv 1(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)$
$G(x)A(x)\equiv 1(\mathrm{m}\mathrm{o}\mathrm{d}{x}^{\frac{n}{2}})$
同理
$$\begin{array}{rl}B(x)-G(x)& \equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^{\frac{n}{2}})\\ (B(x)-G(x){)}^2& \equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)\\ B^2(x)-2B(x)G(x)+G^2(x)& \equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)\\ B(x)-2G(x)+G^2(x)A(x)& \equiv 0(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)\\ B(x)& \equiv 2G(x)-A(x)G^2(x)(\mathrm{m}\mathrm{o}\mathrm{d}{x}^n)\end{array}$$
答案的常数项就是A的常数项的逆元

---

Details
求你开根的次数界要开两倍，每次做完要将$n-2n$的删掉

#include<cstring>
#include<cstdio>
#include<algorithm>
#define mo 998244353
#define REP(i,a,b) for(int i=a,___=b;i<___;++i)
#define fo(i,a,b) for(int i=a,___=b;i<=___;++i)
#define fd(i,a,b) for(int i=a,___=b;i>=___;--i)
#define M 600100
#define n2 499122177
#define ny 332748118
#define L 131071
#define P 524288
#define ll long long
using namespace std;

int G[M],n,N,a[M],b[M],t[M],T[M],w[M],rev[M],fr[M],ifr[M],tri[M],itri[M];

ll sqr(ll a){return a*a%mo;}
ll qpow(ll a,ll i){
	ll r=1;for(;i;i>>=1,a=a*a%mo)if(i&1)r=r*a%mo;return r;
}

void ntt(int *a){
	REP(i,1,N)rev[i]=(rev[i>>1]>>1)+(i&1)*(N>>1),swap(a[i],a[max(rev[i],i)]);
	for(int h=1,m=2;h<N;h=m,m<<=1)REP(i,0,h)for(int j=i,k,T;j<N;j+=m)
		k=j+h,T=1ll*w[P/m*i]*a[k]%mo,a[k]=(a[j]+mo-T)%mo,a[j]=(a[j]+T)%mo;
}
void intt(int *a){
	REP(i,1,N>>1)swap(a[i],a[N-i]);ntt(a);
	ll n=qpow(N,mo-2);REP(i,0,N)a[i]=n*a[i]%mo;
}
void inv(int *A){
	REP(i,0,N<<1)a[i]=T[i]=0;
	a[0]=qpow(A[0],mo-2);int n=N;
	for(N=2;N<=n;){
		REP(i,0,N)T[i]=A[i];
		N<<=1;ntt(T);ntt(a);
		REP(i,0,N)a[i]=(2ll*a[i]+mo-1ll*a[i]*a[i]%mo*T[i]%mo)%mo;
		intt(a);REP(i,N>>1,N)a[i]=0;
	}N>>=1;REP(i,0,N)A[i]=a[i];
}
void sqroot(int *A){
	REP(i,0,N<<1)b[i]=t[i]=0;
	b[0]=1;int n=N;
	for(N=2;N<=n;){
		REP(i,0,N)t[i]=b[i];inv(t);
		REP(i,0,N)T[i]=A[i];REP(i,N,N<<1)T[i]=0;
		N<<=1;ntt(t);ntt(T);ntt(b);
		REP(i,0,N)b[i]=(1ll*T[i]*n2%mo*t[i]+1ll*b[i]*n2)%mo;
		intt(b);REP(i,N>>1,N)b[i]=0;
	}N>>=1;REP(i,0,N)A[i]=b[i];
}

int main(){
	ll W=qpow(3,(mo-1)/P);w[0]=1;REP(i,1,P)w[i]=W*w[i-1]%mo;
	fr[0]=ifr[0]=tri[0]=itri[0]=1;fo(i,1,L)fr[i]=1ll*fr[i-1]*i%mo;
	ifr[L]=qpow(fr[L],mo-2);fd(i,L-1,1)ifr[i]=1ll*ifr[i+1]*(i+1)%mo;
	fo(i,1,L)tri[i]=3ll*tri[i-1]%mo,itri[i]=1ll*ny*itri[i-1]%mo;
	fo(i,0,L)if(i&1)a[i]=sqr(qpow(itri[i>>1],i+1>>1))*ifr[i]%mo;
			 else b[i]=sqr(qpow(itri[i>>1],i>>1))*ifr[i]%mo;
	N=262144;ntt(a);ntt(b);
	REP(i,0,N)G[i]=(1ll*a[i]*a[i]%mo*ny+1ll*b[i]*b[i]+2ll*a[i]*b[i])%mo;
	intt(G);
	REP(i,0,N)if(i&1)G[i]=sqr(qpow(tri[i>>1],i+1>>1))*G[i]%mo;
			 else G[i]=sqr(qpow(tri[i>>1],i>>1))*G[i]%mo;
	sqroot(G);
	scanf("%d",&n);fo(i,1,n)scanf("%d",&n),printf("%d\n",1ll*G[n]*fr[n]%mo);
}