正题

题目链接:https://www.luogu.com.cn/problem/P8338

---

题目大意

给出一个排列$p_i$，定义$a_i^0=i,a_i^k=a_{p_i}^{k-1}$。

对于排列$P$，定义$F(P)$表示最小的一个正整数$k$满足$P^{k+1}=P$。
定义$f(i,j)$，若存在一个$p_i^k=p_j$那么$f(i,j)=0$，否则记$P^{\prime }$表示将$p_i$和$p_j$交换后的排列，$f_{i,j}=F(P^{\prime })$

求
$$\sum _{i=1}^n\sum _{j=1}^{n}f(i,j)$$

答案对$10^9+7$取模

$1\le n\le 5\times 10^5,1\le T\le 5$

---

解题思路

考虑置换环，$i\to p_i$，那么$F(P)$就是所有环的大小的LCM。
然后交换$p_i,p_j$的话就相当于把两个环拼到一起。

注意到所有环的长度和为$n$，那么就证明不同的长度不超过$\sqrt{n}$种，我们可以考虑同一种长度一起处理。

因为最多删去两个环，处理出每个环长质因数幂数最大的前三个。

然后暴力处理就好了，需要记得提前预处理处所有数的质因数分解。

时间复杂度：$O(Tn\log n)$

---

code

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<vector>
#include<cctype>
#include<stack>
#define ll long long
#define mp(x,y) make_pair(x,y)
using namespace std;
const ll N=5e5+10,P=1e9+7;
ll T,n,m,ans,p[N],c[N];
ll inv[N],val[N],fa[N];
bool v[N],vis[N];
pair<ll,ll> a[N];
vector<pair<ll,ll> >q[N];
vector<ll>z[N],pw[N];
stack<ll> cl;
ll read(){ll x=0,f=1;char c=getchar();while(!isdigit(c)){if(c=='-')f=-f;c=getchar();}while(isdigit(c)){x=(x<<1)+(x<<3)+c-'0';c=getchar();}return x*f;
}
ll dfs(ll x){if(v[x])return 0;v[x]=1;return dfs(fa[x])+1;
}
void Ins(ll x){for(ll i=0;i<q[x].size();i++){ll a=q[x][i].first,b=q[x][i].second;z[a].push_back(b);}
}
void Del(ll x){for(ll i=0;i<q[x].size();i++){ll a=q[x][i].first,b=q[x][i].second;if(b==z[a][p[a]]){ans=ans*inv[pw[a][z[a][p[a]]]]%P;if(!p[a])cl.push(a);p[a]++;ans=ans*pw[a][z[a][p[a]]]%P;}}return;
}
void Add(ll x){for(ll i=0;i<q[x].size();i++){ll a=q[x][i].first,b=q[x][i].second;if(b>z[a][p[a]]){ans=ans*inv[pw[a][z[a][p[a]]]]%P;ans=ans*pw[a][b]%P;}}
}
bool cmp(ll x,ll y){return x>y;}
signed main()
{
//	freopen("perm3.in","r",stdin);T=read();inv[0]=inv[1]=1;for(ll i=2;i<N;i++)inv[i]=P-inv[P%i]*(P/i)%P;for(ll i=1;i<N;i++)val[i]=i;for(ll i=2;i<N;i++){if(!vis[i]){pw[i].push_back(1);for(ll j=i;j<N;j=j*i)pw[i].push_back(j);for(ll j=i+i;j<N;j+=i)vis[j]=1;}for(ll j=i;j<N;j+=i){if(val[j]%i==0){ll c=0;while(val[j]%i==0)val[j]/=i,c++;q[j].push_back(mp(i,c));}}}while(T--){m=0;n=read();for(ll i=1;i<=n;i++)z[i].clear(),p[i]=v[i]=c[i]=0;for(ll i=1;i<=n;i++)fa[i]=read();for(ll i=1;i<=n;i++)if(!v[i])c[dfs(i)]++;for(ll i=1;i<=n;i++)if(c[i]){a[++m]=mp(i,c[i]);for(ll j=1;j<=c[i];j++)Ins(i);}ans=1;for(ll i=2;i<=n;i++)if(!vis[i]){z[i].push_back(0);sort(z[i].begin(),z[i].end(),cmp);p[i]=0;ans=ans*pw[i][z[i][0]]%P;}ll pre=ans,sum=0;for(ll i=1;i<=m;i++){for(ll j=1;j<i;j++){Del(a[i].first);Del(a[j].first);Add(a[i].first+a[j].first);(sum+=2ll*ans*a[i].second*a[i].first%P*a[j].second*a[j].first%P)%=P;while(!cl.empty())p[cl.top()]=0,cl.pop();ans=pre;}if(a[i].second>1){Del(a[i].first);Del(a[i].first);Add(a[i].first*2);(sum+=ans%P*a[i].second*a[i].first%P*(a[i].second-1)*a[i].first%P)%=P;while(!cl.empty())p[cl.top()]=0,cl.pop();ans=pre;}}printf("%lld\n",sum);}return 0;
}