这是一篇简洁易懂的良心题解，提供了多种做法。
对于两个数 \(x,y\)，定义 \(x=n^2 \cdot tx,y=m^2 \cdot ty\)。如果 \(x\cdot y\) 为平方数，则说明 \(tx=ty\)。所以我们可以将每个数除去其平方因子，比较相邻两数是否相等即可。
F1：
定义 \(f_{i,j,k}\) 为插入 \(i\) 个数、有 \(j\) 对相邻相等的数、有 \(k\) 对相邻且与 \(a_i\) 相等的数。
定义 \(s_i\) 为在插入第 \(i\) 个数之前有多少个数与 \(a_i\) 相等。
转移分三种情况讨论：
- \(a_i\) 插入后与等于 \(a_i\) 的数相邻：
- \(a_i\) 插入后相邻的两个数相等：
- \(a_i\) 插入后没发生任何变化：
实现时可以把 \(a\) 数组排序，如果 \(a_i \neq a_{i-1}\)，则将 \(f_{i-1,j,k}\) 全部加到 \(f_{i-1,j,0}\) 上，并清 \(0\)。
可以把第一维滚掉。复杂度 \(O(n^3)\)。
code:
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=300+5,mod=1e9+7;
ll n,a[N],f[2][N][N];
void fen(ll &x){
for(ll i=2;i*i<=x;++i){
while(x%(i*i)==0){
x/=i*i;
}
}
}
int main(){
ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
cin>>n;for(int i=1;i<=n;++i){cin>>a[i];fen(a[i]);}
sort(a+1,a+n+1);
int t=0,cnt=0;f[0][0][0]=1;
for(int i=1;i<=n;++i){
if(a[i]!=a[i-1]){
for(int j=0;j<i;++j){
for(int k=1;k<=min(j,cnt);++k){
f[t^1][j][0]=(f[t^1][j][0]+f[t^1][j][k])%mod;
f[t^1][j][k]=0;
}
}
cnt=0;
}
for(int j=0;j<i;++j){
for(int k=0;k<=min(cnt,j);++k){
ll tmp=0;
if(j&&k)tmp=f[t^1][j-1][k-1]*(2*cnt-k+1)%mod;
tmp=(tmp+f[t^1][j][k]*(i-2*cnt+2*k-j)%mod+f[t^1][j+1][k]*(j-k+1)%mod)%mod;
f[t][j][k]=(f[t][j][k]+tmp)%mod;
}
}
++cnt;memset(f[t^1],0,sizeof(f[t^1]));t^=1;
}
cout<<f[t^1][0][0]<<endl;
return 0;
}
F2：
考虑将相同的数作为一个整体，切成多块插入。
定义 \(f_{i,j}\) 为前 \(i\) 组数插入后有 \(j\) 对相同相邻数的方案数。
定义 \(g_{i,j}\) 为 \(i\) 个不同的数分成 \(j\) 组的方案数。
定义 \(s_i\) 为第 \(i\) 组有多少数，\(k\) 为将第 \(i\) 组数分成几部分插入，\(sum_i\) 为 \(\sum\limits_{i=1}^{i}s_i\)，\(p\) 为插入 \(k\) 部分后有几对原来相邻相同的数断开了。
\(g\) 的转移方程如下：
第一项为将 \(i\) 放到块中的方案数，第二项为将 \(i\) 独立成新的一块的方案数。
\(f\) 的转移方程如下：
复杂度 \(O(n^2\sum{s_i})\)，即 \(O(n^3)\)。
code:
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=300+10,mod=1e9+7;
ll n,m,s[N],sum[N],a[N],f[N][N],g[N][N],c[N][N];
void fen(ll &x){
for(ll i=2;i*i<=x;++i){
while(x%(i*i)==0){
x/=i*i;
}
}
}
int main(){
ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
cin>>n;for(int i=1;i<=n;++i){cin>>a[i];fen(a[i]);}
sort(a+1,a+n+1);
for(int i=1;i<=n;++i){
m+=(a[i]!=a[i-1]);
++s[m];
}
for(int i=1;i<=m;++i)sum[i]=sum[i-1]+s[i];
g[0][0]=f[0][0]=c[0][0]=1;
for(int i=1;i<=n+5;++i){
c[i][0]=1;
for(int j=1;j<=i;++j){
g[i][j]=((i-1+j)*g[i-1][j]%mod+j*g[i-1][j-1]%mod)%mod;
c[i][j]=(c[i-1][j-1]+c[i-1][j])%mod;
}
}
for(int i=1;i<=m;++i){
for(int j=0;j<=max(sum[i-1]-1,0LL);++j){
for(int k=0;k<=s[i];++k){
for(int p=0;p<=min(k,j);++p){
f[i][j+s[i]-k-p]=(f[i][j+s[i]-k-p]+c[j][p]*c[sum[i-1]+1-j][k-p]%mod*g[s[i]][k]%mod*f[i-1][j]%mod)%mod;
}
}
}
}
cout<<f[m][0]<<endl;
return 0;
}
F3：
神仙容斥。
我们同样是考虑将相同数字分为一组，对于每一组再分块。
定义 \(s_i\) 为第 \(i\) 组元素个数，\(b_i\) 为分成几块，\(B\) 为 \(\sum\limits_{i=1}^{m}b_i\)。
那么第 \(i\) 组分块排列的方案数为：\(s_i!\cdot\dbinom{s_i-1}{b_i-1}\)。即全排列后插板分块。
那么答案的计算式如下：
因为 \(n,m,s_i\) 是已知的，我们提到前面整理一下：
容易发现，在 \(B\) 确定的情况下，后面那些是可以预处理的。
定义 \(f_{i,j}=\sum{b}[\sum\limits_{w=1}^{i}b_w=j]\left(\prod\limits_{w=1}^{i}\dbinom{s_w-1}{b_w-1}\cdot\frac{1}{b_w!}\right)\)。
转移可以枚举 \(b_i\) 大小，方程如下：
所以最终答案为：
预处理复杂度 \(O(n\sum s_i)\)，即 \(O(n^2)\)。查询复杂度 \(O(n)\)。
总复杂度 \(O(n^2)\)。
code:
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=300+5,mod=1e9+7;
ll n,m,a[N],s[N],fac[N],inv[N],f[N][N];
void fen(ll &x){for(ll i=2;i*i<=x;++i)while(x%(i*i)==0)x/=i*i;}
ll C(int x,int y){
if(y>x||y<0)return 0;
return fac[x]*inv[x-y]%mod*inv[y]%mod;
}
int main(){
ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
cin>>n;for(int i=1;i<=n;++i){cin>>a[i];fen(a[i]);}
sort(a+1,a+n+1);
fac[0]=fac[1]=inv[0]=inv[1]=1;
for(int i=2;i<=n;++i){
fac[i]=fac[i-1]*i%mod;
inv[i]=(mod-mod/i)*inv[mod%i]%mod;
}
for(int i=2;i<=n;++i)inv[i]=inv[i-1]*inv[i]%mod;
for(int i=1;i<=n;++i){
m+=(a[i]!=a[i-1]);
++s[m];
}
ll t=1,sum=0,ans=0;
for(int i=1;i<=m;++i)t=t*fac[s[i]]%mod;
f[0][0]=1;
for(int i=1;i<=m;++i){
sum+=s[i];
for(int j=1;j<=sum;++j){
for(int k=1;k<=min((ll)j-i+1,s[i]);++k){
f[i][j]=(f[i][j]+f[i-1][j-k]*C(s[i]-1,k-1)%mod*inv[k]%mod)%mod;
}
}
}
for(int i=m;i<=n;++i){
ll tmp=t*fac[i]%mod*f[m][i]%mod,op=1;
if((n-i)&1)op=-1;
ans=(ans+op*tmp+mod)%mod;
}
cout<<ans<<endl;
return 0;
}