很早之前看到过这道题，当时连置换是什么都不知道。。。
首先考虑置换ggg中的某一个长为SSS的循环环，学过群论的应该知道，在置换自乘nnn次之后，这个循环会变为gcd(n,S)gcd(n,S)gcd(n,S)个循环。
考虑将ttt个长为SSS的循环拼接成一个，使得若干次自乘之后这个循环断裂为ttt个长为SSS的循环。翻过任何一本群论教材就知道，最后在同一循环中的元素是那些只与初始下标%S\%S%S的值有关，所以这一部分的方案数为St−1(t−1)!S^{t-1}(t-1)!St−1(t−1)!，其实就是考虑一个圆排列。
所有长度为SSS的循环显然可以放在一起处理。
考虑expexpexp的组合意义，显然我们对于每个长度求出拿ttt个该长度的串来拼接成一个的方案数的EGFEGFEGF，设为A(x)A(x)A(x)，则拿出cnt[S]cnt[S]cnt[S]个循环，拼接成若干个的方案数的生成函数为exp(A(x))\exp(A(x))exp(A(x))，复杂度O(nlogn)O(n\logn)O(nlogn)，而且不好卡满。
然而这个玩意可以直接DPDPDP，复杂度O(n43)O(n^{\frac{4}{3}})O(n34)，同样，几乎卡不满
代码（DP）：
#include<bits/stdc++.h>#definelllonglong#definereregister#definegcget_char#definecsconstnamespaceIO{inlinecharget_char(){staticcsintRlen=1<<22|1;staticcharbuf[Rlen],*p1,*p2;return(p1==p2)&&(p2=(p1=buf)+fread(buf,1,Rlen,stdin),p1==p2)?EOF:*p1++;}template<typenameT>inlineTget(){charc;while(!isdigit(c=gc()));Tnum=c^48;while(isdigit(c=gc()))num=(num+(num<<2)<<1)+(c^48);returnnum;}inlineintgetint(){returnget<int>();}}usingnamespaceIO;usingstd::cerr;usingstd::cout;csintmod=998244353;inlineintadd(inta,intb){return(a+=b)>=mod?a-mod:a;}inlineintdec(inta,intb){return(a-=b)<0?a+mod:a;}inlineintmul(inta,intb){llr=(ll)a*b;returnr>=mod?r%mod:r;}inlineintpower(inta,intb,intres=1){for(;b;b>>=1,a=mul(a,a))(b&1)&&(res=mul(res,a));returnres;}inlinevoidInc(int&a,intb){(a+=b)>=mod&&(a-=mod);}inlinevoidDec(int&a,intb){(a-=b)<0&&(a+=mod);}inlineintgcd(inta,intb){#definectz__builtin_ctzintshift=ctz(a|b);for(b>>=ctz(b);a;a-=b)if((a>>=ctz(a))<b)std::swap(a,b);returnb<<=shift;}csintN=1e5+5;intn;intfac[N],ifac[N];inlinevoidinit(){fac[0]=fac[1]=ifac[0]=1;for(intrei=2;i<=n;++i)fac[i]=mul(fac[i-1],i);ifac[n]=power(fac[n],mod-2);for(intrei=n-1;i;--i)ifac[i]=mul(ifac[i+1],i+1);}inlineintC(intn,intm){return(n>=0&&m>=0&&n>=m)?mul(fac[n],mul(ifac[m],ifac[n-m])):0;}std::vector<int>factor;intto[N],cnt[N];boolvis[N];intdp[N],coef[N];signedmain(){#ifdefzxyoifreopen("message.in","r",stdin);#endifn=getint();init();for(intrei=1;i<=n;++i)to[i]=getint();for(intrei=1;i<=n;++i)if(!vis[i]){intlen=0,x=i;while(!vis[x]){vis[x]=true;++len;x=to[x];}++cnt[len];}for(intrei=1;i*i<=n;++i)if(n%i==0){factor.push_back(i);if(i*i!=n)factor.push_back(n/i);}std::sort(factor.begin(),factor.end());intans=1;for(intrei=1;i<=n&&ans;++i)if(cnt[i]){coef[1]=1;for(intrej=2;j<=cnt[i];++j)coef[j]=mul(coef[j-1],mul(i,j-1));memset(dp+1,0,sizeof(int)*cnt[i]);dp[0]=1;for(intrej=1;j<=cnt[i];++j){for(intret:factor){if(t>j)break;if(gcd(n,i*t)!=t)continue;Inc(dp[j],mul(dp[j-t],mul(C(j-1,t-1),coef[t])));}}ans=mul(ans,dp[cnt[i]]);}cout<<ans<<"\n";return0;}
代码exp：
#include<bits/stdc++.h>#definelllonglong#definereregister#definegcget_char#definecsconstnamespaceIO{inlinecharget_char(){staticcsintRlen=1<<22|1;staticcharbuf[Rlen],*p1,*p2;return(p1==p2)&&(p2=(p1=buf)+fread(buf,1,Rlen,stdin),p1==p2)?EOF:*p1++;}template<typenameT>inlineTget(){charc;while(!isdigit(c=gc()));Tnum=c^48;while(isdigit(c=gc()))num=(num+(num<<2)<<1)+(c^48);returnnum;}inlineintgetint(){returnget<int>();}}usingnamespaceIO;usingstd::cerr;usingstd::cout;inlineintgcd(inta,intb){#definectz__builtin_ctzintshift=ctz(a|b);for(b>>=ctz(b);a;a-=b)if((a>>=ctz(a))<b)std::swap(a,b);returnb<<=shift;}csintmod=998244353;inlineintadd(inta,intb){return(a+=b)>=mod?a-mod:a;}inlineintdec(inta,intb){return(a-=b)<0?a+mod:a;}inlineintmul(inta,intb){llr=(ll)a*b;returnr>=mod?r%mod:r;}inlineintpower(inta,intb,intres=1){for(;b;b>>=1,a=mul(a,a))(b&1)&&(res=mul(res,a));returnres;}inlinevoidInc(int&a,intb){(a+=b)>=mod&&(a-=mod);}inlinevoidDec(int&a,intb){(a-=b)<0&&(a+=mod);}csintSIZE=1<<19|1,bit=19;typedefstd::vector<int>Poly;intr[SIZE],*w[bit+1];inlinevoidinit_NTT(){for(intrei=1;i<=bit;++i)w[i]=newint[1<<i-1];intwn=power(3,(mod-1)>>bit);w[bit][0]=1;for(intrei=1;i<(1<<bit-1);++i)w[bit][i]=mul(w[bit][i-1],wn);for(intrei=bit-1;i;--i)for(intrej=0;j<(1<<i-1);++j)w[i][j]=w[i+1][j<<1];}inlinevoidNTT(Poly&A,intlen,inttyp){for(intrei=0;i<len;++i)if(i<r[i])std::swap(A[i],A[r[i]]);for(intrei=1,t=1;i<len;i<<=1,++t)for(intrej=0;j<len;j+=i<<1)for(intrek=0;k<i;++k){intx=A[j+k],y=mul(A[j+k+i],w[t][k]);A[j+k]=add(x,y),A[j+k+i]=dec(x,y);}if(typ==-1){std::reverse(A.begin()+1,A.begin()+len);for(intrei=0,inv=power(len,mod-2);i<len;++i)A[i]=mul(A[i],inv);}}inlinevoidinit_rev(intlen){for(intrei=0;i<len;++i)r[i]=r[i>>1]>>1|((i&1)?len>>1:0);}intfac[SIZE],ifac[SIZE],inv[SIZE];inlinevoidinit_inv(){fac[0]=fac[1]=ifac[0]=ifac[1]=inv[0]=inv[1]=1;for(intrei=2;i<SIZE;++i){fac[i]=mul(fac[i-1],i);inv[i]=mul(inv[mod%i],mod-mod/i);ifac[i]=mul(ifac[i-1],inv[i]);}}inlineintC(intn,intm){return(n>=0&&m>=0&&n>=m)?mul(fac[n],mul(ifac[m],ifac[n-m])):0;}inlinePolyoperator*(Polya,Polyb){intdeg=a.size()+b.size()-1,l=1;while(l<deg)l<<=1;init_rev(l);a.resize(l);NTT(a,l,1);b.resize(l);NTT(b,l,1);for(intrei=0;i<l;++i)a[i]=mul(a[i],b[i]);NTT(a,l,-1);a.resize(deg);returna;}inlinePolyDeriv(Polya){for(intrei=0;i+1<a.size();++i)a[i]=mul(a[i+1],i+1);a.pop_back();returna;}inlinePolyInteg(Polya){a.push_back(0);for(intrei=a.size()-1;i;--i)a[i]=mul(a[i-1],inv[i]);a[0]=0;returna;}inlinePolyInv(csPoly&a,intlim){Polyc,b(1,power(a[0],mod-2));for(intrel=4;(l>>2)<lim;l<<=1){init_rev(l);c=a,c.resize(l>>1);c.resize(l),NTT(c,l,1);b.resize(l),NTT(b,l,1);for(intrei=0;i<l;++i)b[i]=mul(b[i],dec(2,mul(b[i],c[i])));NTT(b,l,-1);b.resize(l>>1);}b.resize(lim);returnb;}inlinePolyLn(Polya,intlim){a=Integ(Deriv(a)*Inv(a,lim));a.resize(lim);returna;}inlinePolyExp(csPoly&a){Polyc,b(1,1);intn=a.size();for(intrei=2;(i>>1)<n;i<<=1){c=Ln(b,i);for(intrej=0;j<i;++j)c[j]=dec(j<n?a[j]:0,c[j]);c[0]=add(c[0],1);b=b*c;b.resize(i);}b.resize(n);returnb;}csintN=1e5+5;intn;std::vector<int>factor;intto[N],cnt[N];boolvis[N];Polyf,g;signedmain(){init_NTT();init_inv();#ifdefzxyoifreopen("message.in","r",stdin);#endifn=getint();for(intrei=1;i<=n;++i)to[i]=getint();for(intrei=1;i<=n;++i)if(!vis[i]){intlen=0,x=i;while(!vis[x]){vis[x]=true;++len;x=to[x];}++cnt[len];}for(intrei=1;i*i<=n;++i)if(n%i==0){factor.push_back(i);if(i*i!=n)factor.push_back(n/i);}std::sort(factor.begin(),factor.end());intans=1;for(intrei=1;i<=n&&ans;++i)if(cnt[i]){f=Poly(cnt[i]+1,0);for(intrej:factor){if(j>cnt[i])break;if(gcd(j*i,n)!=j)continue;f[j]=power(i,j-1,inv[j]);}g=Exp(f);ans=mul(ans,mul(g[cnt[i]],fac[cnt[i]]));}cout<<ans<<"\n";return0;}
