函树(hs) - 莫比乌斯反演 - 虚树

题目大意：给一颗树，求${\sum }_{x,y}\phi (x\times y)\text{dist}(x,y)$。$n\le 10^5$
题解：
$$\begin{gathered} \sum _{x,y}\phi (x\times y)\text{dist}(x,y) \\ =\sum _{x,y}\frac{\phi (x)\phi (y)\gcd (x,y)}{\phi (\gcd (x,y))}\text{dist}(x,y) \\ =\sum _{d=1}^n\frac{d}{\phi (d)}\sum _{d|x,d|y}\phi (x)\phi (y)\text{dist}(x,y)\left[\gcd (x,y)=d\right] \\ =\sum _{d=1}^n\frac{d}{\phi (d)}\sum _{d|x,d|y}\phi (x)\phi (y)\text{dist}(x,y)\sum _{e|\frac{x}{d},e|\frac{y}{d}}\mu (e) \\ =\sum _{T=1}^n\text{F}(T)\sum _{d|T}\frac{d}{\phi (d)}\mu \left(\frac{T}{d}\right) \end{gathered}$$
其中$\text{F}(T)$表示把所有两端点编号是$T$的倍数的路径拿出来，求端点的$\phi$的乘积再求和。这个直接虚树即可。实现优秀可以$O\left(n\lg n\right)$

#include<bits/stdc++.h>
#define rep(i,a,b) for(int i=a;i<=b;i++)
#define Rep(i,v) rep(i,0,(int)v.size()-1)
#define lint long long
#define mod 998244353
#define ull unsigned lint
#define db long double
#define pb push_back
#define mp make_pair
#define fir first
#define sec second
#define debug(x) cerr<<#x<<"="<<x
#define sp <<" "
#define ln <<endl
using namespace std;
typedef pair<int,int> pii;
typedef set<int>::iterator sit;
namespace INPUT_SPACE{
    const int BS=(1<<24)+5;char Buffer[BS],*HD,*TL;inline int gc() { if(HD==TL) TL=(HD=Buffer)+fread(Buffer,1,BS,stdin);return (HD==TL)?EOF:*HD++; }
    inline int inn() { int x,ch;while((ch=gc())<'0'||ch>'9');x=ch^'0';while((ch=gc())>='0'&&ch<='9') x=(x<<1)+(x<<3)+(ch^'0');return x; }
}using INPUT_SPACE::inn;
const int N=100010;
struct edges{
    int to,pre;
}e[N<<1],dvs[1166750+100];int h[N],etop,dh[N],dtop,*lst[N],lstc[N],d[N],tms[N];
inline int add_edge(int u,int v) { return e[++etop].to=v,e[etop].pre=h[u],h[u]=etop; }
inline int add_dvs(int u,int v) { return dvs[++dtop].to=v,dvs[dtop].pre=dh[u],dh[u]=dtop; }
inline int fast_pow(int x,int k,int ans=1) { for(;k;k>>=1,x=(lint)x*x%mod) if(k&1) ans=(lint)ans*x%mod;return ans; }
namespace PRELUDE_space{
    int mu[N],phi[N],val[N],p[N],np[N],fac[N],facinv[N],inv[N];
    inline int prelude(int n)
    {
        phi[1]=mu[1]=1;
        for(int i=2,c=0;i<=n;i++)
        {
            if(!np[i]) p[++c]=i,phi[i]=i-1,mu[i]=-1;
            rep(j,1,c&&p[j]<=n/i)
            {
                int x=i*p[j];np[x]=1;
                if(i%p[j]==0) { phi[x]=phi[i]*p[j],mu[x]=0;break; }
                else phi[x]=phi[i]*(p[j]-1),mu[x]=-mu[i];
            }
        }
        rep(i,fac[0]=1,n) fac[i]=(lint)fac[i-1]*i%mod;
        facinv[n]=fast_pow(fac[n],mod-2);
        for(int i=n-1;i>=0;i--) facinv[i]=(i+1ll)*facinv[i+1]%mod;
        rep(i,1,n) inv[i]=(lint)facinv[i]*fac[i-1]%mod;
        rep(d,1,n) rep(e,1,n/d) if(mu[e])
        {
            int t=d*e;
            if(mu[e]>0) val[t]=(val[t]+(lint)d*inv[phi[d]])%mod;
            else val[t]=(val[t]-(lint)d*inv[phi[d]])%mod,(val[t]<0?val[t]+=mod:0);
        }
        return 0;
    }
}using PRELUDE_space::val;using PRELUDE_space::phi;
namespace LCA_space{
    const int LOG=20;
    int in[N],fir[N],f[N<<1][LOG],Log[N<<1],dfc,cnt;
    int dfs(int x,int fa=0)
    {
        tms[in[x]=++dfc]=x,d[x]=d[fa]+1,f[fir[x]=++cnt][0]=x;
        for(int i=h[x],y;i;i=e[i].pre)
            if((y=e[i].to)^fa) dfs(y,x),f[++cnt][0]=x;
        return 0;
    }
    inline int Mymin(int x,int y) { return d[x]<d[y]?x:y; }
    inline int build_ST()
    {
        rep(i,2,cnt) Log[i]=Log[i>>1]+1;
        for(int j=1;(1<<j)<=cnt;j++) rep(i,1,cnt-(1<<j)+1)
            f[i][j]=Mymin(f[i][j-1],f[i+(1<<(j-1))][j-1]);
        return 0;
    }
    inline int LCA(int x,int y)
    {
        x=fir[x],y=fir[y];if(x>y) swap(x,y);int k=Log[y-x+1];
        return Mymin(f[x][k],f[y-(1<<k)+1][k]);
    }
}using LCA_space::dfs;using LCA_space::LCA;
namespace SOLVE_space{
    int p[N],dp[N],imp[N],ans,s[N];vector<int> g[N];
    inline int getans(int x)
    {
        p[x]=dp[x]=0;if(imp[x]) p[x]=phi[x];
        Rep(i,g[x])
        {
            int y=g[x][i];getans(y);
            dp[y]=(dp[y]+(lint)p[y]*(d[y]-d[x]))%mod,
            ans=(ans+(lint)p[x]*dp[y]+(lint)p[y]*dp[x])%mod,
            dp[x]+=dp[y],(dp[x]>=mod?dp[x]-=mod:0),
            p[x]+=p[y],(p[x]>=mod?p[x]-=mod:0);
        }
        return 0;
    }
    inline int solve(int *lst,int n,int qwq)
    {
        int t;s[t=1]=lst[1];
        rep(i,1,n) g[lst[i]].clear();
        rep(i,2,n)
        {
            int x=lst[i],y=s[t],z=s[t-1],c;
            while(d[c=LCA(x,y)]<=d[z]) g[z].pb(y),y=s[--t],z=s[t-1];
            if(c^y) g[c].clear(),g[c].pb(y),s[t]=c;s[++t]=x;
        }
        while(t>1) g[s[t-1]].pb(s[t]),t--;
        rep(i,1,n) imp[lst[i]]=1;imp[1]=qwq;
        ans=0,getans(lst[1]);
        rep(i,1,n) imp[lst[i]]=0;
        return ans;
    }
}using SOLVE_space::solve;
int main()
{
    int n=inn(),x,y,ans=0;PRELUDE_space::prelude(n);
    rep(i,1,n-1) x=inn(),y=inn(),add_edge(x,y),add_edge(y,x);
    dfs(1),LCA_space::build_ST();
    rep(i,1,n) rep(j,1,n/i) add_dvs(i*j,i);
    rep(i,1,n) lst[i]=new int[n/i+5];
    rep(i,2,n) lst[i][lstc[i]=1]=1;
    rep(t,1,n)
    {
        int x=tms[t];
        for(int i=dh[x],y;i;i=dvs[i].pre)
            y=dvs[i].to,lst[y][++lstc[y]]=x;
    }
    rep(i,1,n) ans=(ans+2ll*solve(lst[i],lstc[i],i==1)*val[i])%mod;
    return !printf("%d\n",ans);
}