zoukankan      html  css  js  c++  java
  • YbtOJ#723欧拉之树【莫比乌斯反演,虚树】

    正题

    题目链接:http://www.ybtoj.com.cn/contest/121/problem/2


    题目大意

    给出\(n\)个点的一棵树,每个点有一个权值\(a_i\),求

    \[\sum_{i=1}^n\sum_{j=1}^ndis(i,j)\times \varphi(a_i\times a_j) \]

    \(2\leq n\leq 2\times 10^5\)\(a\)恰好是一个排列。


    解题思路

    一个十分显然的结论就是\(\varphi(x\times y)=\varphi(x)\times \varphi(y)\times \frac{gcd(x,y)}{\varphi(gcd(x,y))}\)。(相同的质因子只保留一个数\(p-1\)的就好了)

    然后顺便把点编号换一下使得\(a_i=i\)再枚举约数就是

    \[\sum_{d=1}^n\frac{\varphi(d)}{d}\sum_{i=1}^n\sum_{j=1}^{n}dis(i,j) \varphi(i)\varphi(j)\times [gcd(i,j)=d] \]

    然后就可以莫反了,定义

    \[g_d=\sum_{d|i}^n\sum_{d|j}^ndis(i,j)\varphi(i)\varphi(j) \]

    \[g_d=\sum_{d|i}^n\sum_{d|j}^n(dep_i+dep_j-2dep_{lca(i,j)})\varphi(i)\varphi(j) \]

    \[g_d=2\sum_{d|i}^ndep_{i}\varphi(i)\sum_{d|j}^n\varphi(j)-2\sum_{k=1}^n\sum_{i=1}^n\sum_{j=1}^n[lca(i,j)=k]dep_{k}\varphi(i)\varphi(j) \]

    把所有\(d\)倍的点加入虚树,然后用树形\(dp\)计算后面那个东西,前面那个可以直接算。

    然后答案就是

    \[\sum_{d=1}^n\frac{\varphi(d)}{d}\sum_{d|i}g_i \]

    时间复杂度\(O(n\log^2 n)\),有点卡常。


    code

    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #include<cctype>
    #define ll long long
    #pragma GCC optimize(2)
    %:pragma GCC optimize(3)
    %:pragma GCC optimize("Ofast")
    %:pragma GCC optimize("inline")
    using namespace std;
    const int N=4e5+10,T=20,P=1e9+7;
    int read() {
    	int 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-48,c=getchar();
    	return x*f;
    }
    struct node{
    	int to,next;
    }a[N<<1];
    int n,m,tot,top,dfc,ls[N],rfl[N];
    int mu[N],phi[N],pri[N],prn,s[N],p[N];
    int stn,lg[N],wz[N],rfn[N],dep[N],f[N][T];
    ll S[N],dp[N],g[N],ans;bool v[N],mark[N];
    ll power(ll x,ll b){
    	ll ans=1;
    	while(b){
    		if(b&1)ans=ans*x%P;
    		x=x*x%P;b>>=1;
    	}
    	return ans;
    }
    void prime(){
    	mu[1]=phi[1]=1;
    	for(int i=2;i<=n;i++){
    		if(!v[i])pri[++prn]=i,phi[i]=i-1,mu[i]=-1;
    		for(int j=1;j<=prn&&i*pri[j]<=n;j++){
    			v[i*pri[j]]=1;
    			if(i%pri[j]==0){
    				phi[i*pri[j]]=phi[i]*pri[j];
    				break;
    			}
    			phi[i*pri[j]]=phi[i]*(pri[j]-1);
    			mu[i*pri[j]]=-mu[i];
    		}
    	}
    	return;
    }
    void addl(int x,int y){
    	a[++tot].to=y;
    	a[tot].next=ls[x];
    	ls[x]=tot;return;
    }
    bool cmp(int x,int y)
    {return rfn[x]<rfn[y];}
    void dfs(int x,int fa){
    	dep[x]=dep[fa]+1;rfn[x]=++dfc;
    	f[++stn][0]=x;wz[x]=stn;
    	for(int i=ls[x];i;i=a[i].next){
    		int y=a[i].to;
    		if(y==fa)continue;
    		dfs(y,x);f[++stn][0]=x;
    	}
    	return;
    }
    int LCA(int l,int r){
    	l=wz[l];r=wz[r];
    	if(l>r)swap(l,r);
    	int z=lg[r-l+1],x=f[l][z],y=f[r-(1<<z)+1][z];
    	return dep[x]<dep[y]?x:y;
    }
    void Ins(int x){
    	if(!top){s[++top]=x;return;}
    	int lca=LCA(x,s[top]);
    	while(top>1&&dep[s[top-1]]>dep[lca])
    		addl(s[top-1],s[top]),top--;
    	if(dep[s[top]]>dep[lca])addl(lca,s[top]),top--;
    	if((!top)||s[top]!=lca)s[++top]=lca;
    	s[++top]=x;return;
    }
    void calc(int x,ll &ans){
    	if(mark[x])S[x]=phi[x],dp[x]=1ll*phi[x]*phi[x]%P;
    	else S[x]=dp[x]=0;
    	for(int i=ls[x];i;i=a[i].next){
    		int y=a[i].to;calc(y,ans);
    		(dp[x]+=S[x]*S[y]*2ll%P)%=P;
    		S[x]=(S[x]+S[y])%P;
    	}
    	(ans+=P-1ll*dp[x]*dep[x]%P)%=P;
    	ls[x]=mark[x]=0;return;
    }
    signed main()
    {
    	freopen("sm.in","r",stdin);
    	freopen("sm.out","w",stdout);
    	n=read();prime();
    	for(int i=1;i<=n;i++){
    		int x=read();
    		rfl[i]=x;
    	}
    	for(int i=1;i<n;i++){
    		int x=read(),y=read();
    		x=rfl[x];y=rfl[y];
    		addl(x,y);addl(y,x);
    	}
    	dfs(1,1);
    	for(int j=1;(1<<j)<=stn;j++)
    		for(int i=1;i+(1<<j)-1<=stn;i++){
    			int x=f[i][j-1],y=f[i+(1<<j-1)][j-1];
    			f[i][j]=(dep[x]<dep[y])?x:y;
    		}
    	for(int i=2;i<=stn;i++)lg[i]=lg[i>>1]+1;	
    	memset(ls,0,sizeof(ls));
    	for(int k=1;k<=n;k++){
    		m=top=tot=0;ll sum=0;
    		for(int i=k;i<=n;i+=k)
    			p[++m]=i,sum+=phi[i];
    		sort(p+1,p+1+m,cmp);sum%=P;
    		if(p[1]!=1)s[++top]=1;
    		for(int i=1;i<=m;i++){
    			Ins(p[i]);mark[p[i]]=1;
    			(g[k]+=1ll*phi[p[i]]*dep[p[i]]%P*sum%P)%=P;
    		}
    		while(top>1)addl(s[top-1],s[top]),top--;
    		calc(1,g[k]);g[k]=g[k]*2ll%P;
    	}
    	for(int i=1;i<=n;i++){
    		ll tmp=0;
    		for(int j=i;j<=n;j+=i)
    			(tmp+=mu[j/i]*g[j]%P)%=P;
    		(ans+=tmp*i%P*power(phi[i],P-2)%P)%=P;
    	}
    	printf("%d\n",(ans+P)%P*power(1ll*n*(n-1)%P,P-2)%P);
    	return 0;
    }
    
  • 相关阅读:
    Linux中的samba服务和ftp服务
    Linux中的rpm和yum软件管理
    Linux网络和进程管理
    在jupyter notebook 中编辑公式
    CEO的作用
    如何学习
    敏捷宣言遵循的原则
    敏捷宣言
    python模糊匹配之fuzzywuzzy
    拥抱变革(More Fearless Change)
  • 原文地址:https://www.cnblogs.com/QuantAsk/p/14432453.html
Copyright © 2011-2022 走看看