设$f(x)=sum_{x|d}p(d)$。
则$ans=sum_{i=1}^nsum_{j=1}^nsum_{k=1}^nmu(i)mu(j)mu(k)f(lcm(i,j))f(lcm(i,k))f(lcm(j,k))$。
转化成图论模型，$i$到$j$有边的条件是$mu(i) eq0,mu(j) eq0,lcm(i,j)leq n$。
枚举square-free的$gcd$，再枚举square-free的$lcm$，然后枚举$frac{lcm}{gcd}$的因子$a$，那么可以得到一对满足条件的数对$(a imesgcd,frac{lcm}{a})$。
这样可以不重不漏地枚举出所有边，然后三元环计数即可。
时间复杂度$O(nlog nsqrt{nlog n})$。
#include<cstdio>
#include<algorithm>
using namespace std;
const int N=100005,M=1166760,E=760745;
int n,i,j,k,x,y,z,P[N],f[N],vis[N],mu[N],p[N],tot,ans,A,B;
int g[N],v[M],nxt[M],ed,d[N];
inline void read(int&a){char c;while(!(((c=getchar())>='0')&&(c<='9')));a=c-'0';while(((c=getchar())>='0')&&(c<='9'))(a*=10)+=c-'0';}
inline void add(int x,int y){v[++ed]=y;nxt[ed]=g[x];g[x]=ed;}
namespace Triple{
int g[N],v[E],w[E],nxt[E],ed,st[N],en[N],m,q[M],val[M],tmp[N];
inline void add(int x,int y,int z){
if(d[x]>d[y])swap(x,y);
v[++ed]=y;w[ed]=z;nxt[ed]=g[x];g[x]=ed;
}
void solve(){
for(i=1;i<=n;i++)if(g[i]){
st[i]=m+1;
for(j=g[i];j;j=nxt[j])tmp[q[++m]=v[j]]=w[j];
en[i]=m;
sort(q+st[i],q+m+1);
for(j=st[i];j<=m;j++)val[j]=tmp[q[j]];
}
for(i=1;i<=n;i++)for(j=g[i];j;j=nxt[j]){
k=v[j],x=st[i],y=st[k];
while(x<=en[i]&&y<=en[k])if(q[x]==q[y]){
z=q[x];
A+=mu[i]*mu[k]*mu[z]*w[j]*val[x]*val[y];
x++,y++;
}else q[x]<q[y]?x++:y++;
}
}
}
int main(){
read(n);
for(i=1;i<=n;i++)read(P[i]);
for(i=1;i<=n;i++)for(j=i;j<=n;j+=i)f[i]+=P[j],add(j,i);
for(mu[1]=1,i=2;i<=n;i++){
if(!vis[i])p[tot++]=i,mu[i]=-1;
for(j=0;i*p[j]<=n&&j<tot;j++){
vis[i*p[j]]=1;
if(i%p[j])mu[i*p[j]]=-mu[i];else break;
}
}
for(i=1;i<=n;i++)ans+=mu[i]*f[i]*f[i]*f[i];
for(i=1;i<=n;i++)if(mu[i])for(j=i;j<=n;j+=i)if(mu[j]&&mu[j/i])for(k=g[j/i];k;k=nxt[k]){
x=i*v[k],y=j/v[k];
if(x>=y)continue;
d[x]++,d[y]++;
B+=(mu[x]*f[y]+mu[y]*f[x])*f[j]*f[j];
}
for(i=1;i<=n;i++)if(mu[i])for(j=i;j<=n;j+=i)if(mu[j]&&mu[j/i])for(k=g[j/i];k;k=nxt[k]){
x=i*v[k],y=j/v[k];
if(x>=y)continue;
Triple::add(x,y,f[j]);
}
Triple::solve();
return printf("%d",(ans+A*6+B*3)&((1<<30)-1)),0;
}