经过基本观察，可得当点对 \((i,j)\) 合法时，有 \(a_i|b_j,a_j|b_i\)，其中 \(a_i=i/gcd(p_i,i),b_i=p_i/gcd(p_i,i)\)，证明显然。
如何维护？考虑开 \(mp_{x,y}\) 表示 \(x=a_i\),\(y|b_i\) 的个数。对于点 \(i\) 点对个数即为 \(\sum\limits_{d|b_i} mp_{d,a_i}\)
时间复杂度为 \(O(nlog^2n)\)，空间复杂度为 \(O(nlogn)\)，可通过 \(G1\)。
如何优化？考虑枚举 \(x\)，那么只需一个 \(O(n)\) 的数组维护。令 \(cnt_{(x),y}\) 为 \(a_i=x,y|b_i\) 的点对数，统计答案时枚举 \(x\) 的倍数 \(y\)，
贡献点对个数为 \(\sum\limits_{b_i=y}cnt_{(x),a_i}\)，统计完 \(x\) 后要清空数组。
这样统计会统计到 $(i,j),i\ge j $ 的点对，当 \(a_i|b_i\) 时会被统计，需减掉，根据对称性将最终答案数以 2 即可。
因为每个点只会在 \(p_i\) 的因子处贡献，且 \(\{p_i\}\) 为排列，所以时间复杂度为 \(O(nlogn)\)
//G1
void Main(){
n=rd,ans=0;
for(int i=1;i<=n;i++)
p[i]=rd;
for(int i=1;i<=n;i++){
int g=__gcd(p[i],i);
a[i]=i/g,b[i]=p[i]/g;
// a[i]|b[j],a[j]|b[i]
//统计 ans
for(int j:fac[b[i]])
ans+=mp[j][a[i]];
for(int j:fac[b[i]])
mp[a[i]][j]++;
}
for(int i=1;i<=n;i++)
for(int j:fac[b[i]])
mp[a[i]][j]--;
cout<<ans<<endl;
}
signed main(){
for(int i=1;i<=N-10;i++)
for(int j=i;j<=N-10;j+=i)
fac[j].pb(i);
int T=rd;
while(T--) Main();
}
//G2
void Main(){
n=rd,ans=res=0;
for(int i=1;i<=n;i++)
p[i]=rd;
for(int i=1;i<=n;i++){
int g=__gcd(p[i],i);
a[i]=i/g,b[i]=p[i]/g;
if(b[i]%a[i]==0) res++;
vec1[a[i]].pb(b[i]);
vec2[b[i]].pb(a[i]);
}
for(int x=1;x<=n;x++){
for(int y:vec1[x]) for(int z:fac[y]) cnt[z]++;
for(int y=x;y<=n;y+=x) for(int z:vec2[y]) ans+=cnt[z];
for(int y:vec1[x]) for(int z:fac[y]) cnt[z]--;
}
printf("%lld\n",(ans-res)/2);
for(int i=1;i<=n;i++) vec1[i].clear(),vec2[i].clear();
}