LOJ6271 「长乐集训 2017 Day10」生成树求和 加强版

首先我们可以想到这道题目应该是要用矩阵树定理来解决，但是他需要求的是所有生成树的边权按照三进制不进位加法相加然后再用十进制加法加在一起的和。比较显然的是，我们可以将边权三进制拆分之后每一位单独搞。同时可以利用生成函数的思想，将每一条边的权值设为 \(x^0,x^1,x^2\) 三者中的一种，最后我们求出矩阵树定理后每一项之前的系数即可（注，这里是循环卷积）。

于是我们便大胆地想到，如果我们将边权都变成一个多项式，能否就可以直接套用矩阵树定理做了呢？貌似是可行的，我们可以直接手搓多项式逆。

\[(a_0x^0+a_1x^1+a_2x^2)(b_0x^0+b_1x^1+b_2x^2)=1\\ b_0=\frac{a_0^2-a_1a_2}{a_0^3+a_1^3+a_2^3-3a_0a_1a_2}\\ b_1=\frac{a_2^2-a_0a_1}{a_0^3+a_1^3+a_2^3-3a_0a_1a_2}\\ b_2=\frac{a_1^2-a_0a_2}{a_0^3+a_1^3+a_2^3-3a_0a_1a_2}\\ \]

然后发现存在两种情况可能没有逆，\(a_0+a_1+a_2=0\) 或 \(a_0=a_1=a_2\) 两种情况，但是手摸样例可以发现，如果每次都将没有逆的视作高斯消元中的 \(0\) ，那么最后只有 \((n-1,n-1)\) 的位置才会没有逆，但是这个位置已经不需要存在逆了，所以直接计算即可。

#include<bits/stdc++.h>
using namespace std;
const int N=1e2+5;
const int MOD=1e9+7;
int ADD(int x,int y){return x+y>=MOD?x+y-MOD:x+y;}
int SUB(int x,int y){return x-y<0?x-y+MOD:x-y;}
int TIME(int x,int y){return (int)(1ll*x*y%MOD);}
int ksm(int x,int k=MOD-2){int res=1;for(;k;k>>=1,x=TIME(x,x))if(k&1)res=TIME(res,x);return res;}
int n,m,res=0;
struct Edge{int u,v,w;}a[N*N];
struct Polynomial{int f[3];int &operator [] (int x){return f[x];}};
Polynomial operator * (Polynomial a,Polynomial b){
	Polynomial res;
	res[0]=ADD(TIME(a[0],b[0]),ADD(TIME(a[1],b[2]),TIME(a[2],b[1])));
	res[1]=ADD(TIME(a[0],b[1]),ADD(TIME(a[1],b[0]),TIME(a[2],b[2])));
	res[2]=ADD(TIME(a[0],b[2]),ADD(TIME(a[1],b[1]),TIME(a[2],b[0])));
	return res;
}
Polynomial operator + (Polynomial a,Polynomial b){
	return (Polynomial){ADD(a[0],b[0]),ADD(a[1],b[1]),ADD(a[2],b[2])};
}
Polynomial operator - (Polynomial a,Polynomial b){
	return (Polynomial){SUB(a[0],b[0]),SUB(a[1],b[1]),SUB(a[2],b[2])};
}
Polynomial inv(Polynomial a){
	Polynomial res;int tmp=0;
	tmp=ADD(tmp,TIME(a[0],TIME(a[0],a[0])));
	tmp=ADD(tmp,TIME(a[1],TIME(a[1],a[1])));
	tmp=ADD(tmp,TIME(a[2],TIME(a[2],a[2])));
	tmp=SUB(tmp,TIME(3,TIME(a[0],TIME(a[1],a[2]))));
	tmp=ksm(tmp);
	res[0]=TIME(tmp,SUB(TIME(a[0],a[0]),TIME(a[1],a[2])));
	res[1]=TIME(tmp,SUB(TIME(a[2],a[2]),TIME(a[0],a[1])));
	res[2]=TIME(tmp,SUB(TIME(a[1],a[1]),TIME(a[0],a[2])));
	return res;
}
Polynomial operator / (Polynomial a,Polynomial b){return a*inv(b);}
Polynomial f[N][N];
void print(){
	printf("---------------\n");
	for(int i=1;i<n;++i){
		for(int j=1;j<n;++j){
			printf("[%d,%d,%d] ",f[i][j][0],f[i][j][1],f[i][j][2]);
		}
		printf("\n");
	}
}
int cal(){
	Polynomial res={1,0,0};
	for(int i=1;i<n-1;++i){
		int tmp=i-1;
		for(int j=i;j<n;++j){
			if(ADD(f[j][i][0],ADD(f[j][i][1],f[j][i][2]))==0) continue;
			if(f[j][i][0]==f[j][i][1]&&f[j][i][1]==f[j][i][2]) continue;
			tmp=j;break;
		}
		if(tmp<i) return 0;
		if(tmp!=i){
			res[0]=SUB(0,res[0]),res[1]=SUB(0,res[1]),res[2]=SUB(0,res[2]);
			for(int j=1;j<n;++j) swap(f[i][j],f[tmp][j]);
		}
		res=res*f[i][i];
		for(int j=i+1;j<n;++j){
			Polynomial tmp=f[j][i]/f[i][i];
			for(int k=1;k<n;++k) f[j][k]=f[j][k]-f[i][k]*tmp;
		}
	}
	return res=res*f[n-1][n-1],ADD(res[1],TIME(res[2],2));
}
int main(){
	freopen("sum.in","r",stdin);
	freopen("sum.out","w",stdout);
	cin>>n>>m;
	for(int i=1;i<=m;++i)
		scanf("%d%d%d",&a[i].u,&a[i].v,&a[i].w);
	for(int T=9,tmp=1;T--;tmp*=3){
		for(int i=1;i<=n;++i){
			for(int j=1;j<=n;++j)
			f[i][j]=(Polynomial){0,0,0};
		}
		for(int i=1;i<=m;++i){
			int TMP=a[i].w/tmp%3;
			f[a[i].u][a[i].u][TMP]=ADD(f[a[i].u][a[i].u][TMP],1);
			f[a[i].v][a[i].v][TMP]=ADD(f[a[i].v][a[i].v][TMP],1);
			f[a[i].u][a[i].v][TMP]=SUB(f[a[i].u][a[i].v][TMP],1);
			f[a[i].v][a[i].u][TMP]=SUB(f[a[i].v][a[i].u][TMP],1);
		}
		int TMP=cal();res=ADD(res,TIME(TMP,tmp));
	}
	return printf("%d\n",res),0;
}