【LOJ6271】「长乐集训 2017 Day10」生成树求和 加强版（循环卷积）（矩阵树）

传送门

---

题解：

拆位是肯定要拆位的。

拆位之后需要求的是这一位上三进制不进位加法结果为$0,1,2$的分别有多少。

考虑矩阵树，由于矩阵树求的是边权乘积之和，我们将三进制不进位加法转换成长度为3的循环卷积，求出${\omega }^0,{\omega }^1,{\omega }^2$的点值之后直接IFWT回来，其中$\omega$指三次单位根。

由于${\omega }^2=-\omega -1$，我们可以只用$x+y\times \omega$的复数来表示我们需要的数了。

由于高斯消元需要求逆元，对于$(a+b\omega {)}^{-1}=\frac{1}{a+b\omega }$考虑分母有理化，利用${\omega }^2+\omega =-1$，$\frac{1}{a+b\omega }=\frac{a+b{\omega }^2}{(a+b\omega )(a+b{\omega }^2)}=\frac{(a-b)+(-b)\omega }{a^a+b^2-ab}$

---

代码：

#include<bits/stdc++.h>
#define ll long long
#define re register
#define gc get_char
#define cs const

namespace IO{
	inline char get_char(){
		static cs int Rlen=1<<22|1;
		static char buf[Rlen],*p1,*p2;
		return (p1==p2)&&(p2=(p1=buf)+fread(buf,1,Rlen,stdin),p1==p2)?EOF:*p1++;
	}
	
	template<typename T>
	inline T get(){
		char c;T num;
		while(!isdigit(c=gc()));num=c^48;
		while(isdigit(c=gc()))num=(num+(num<<2)<<1)+(c^48);
		return num;
	}
	inline int gi(){return get<int>();}
}
using namespace IO;

using std::cerr;
using std::cout;

cs int mod=1e9+7,inv3=(mod+1)/3;
inline int add(int a,int b){a+=b-mod;return a+(a>>31&mod);}
inline int dec(int a,int b){a-=b;return a+(a>>31&mod);}
inline int mul(int a,int b){ll r=(ll)a*b;return r>=mod?r%mod:r;}
inline int power(int a,int b,int res=1){
	for(;b;b>>=1,a=mul(a,a))(b&1)&&(res=mul(res,a));
	return res;
}
inline void Inc(int &a,int b){a+=b-mod;a+=a>>31&mod;}
inline void Dec(int &a,int b){a-=b;a+=a>>31&mod;}
inline void Mul(int &a,int b){a=mul(a,b);}
inline void ex_gcd(int a,int b,ll &x,ll &y){
	if(!b){x=1;y=0;return ;}ex_gcd(b,a%b,y,x);y-=a/b*x;
}
inline int inv(int a){
	assert(a);ll x,y;
	ex_gcd(a,mod,x,y);
	return (x%mod+mod)%mod;
}

struct cp{
	int x,y;//x + y\omega_3
	cp(){}
	cp(int _x,int _y=0):x(_x),y(_y){}
	operator bool()cs{return x||y;}
	friend cp operator+(cs cp &a,cs cp &b){return cp(add(a.x,b.x),add(a.y,b.y));}
	friend cp operator-(cs cp &a,cs cp &b){return cp(dec(a.x,b.x),dec(a.y,b.y));}
	friend cp operator*(cs cp &a,cs cp &b){
		ll t=mul(a.y,b.y);return cp(dec(mul(a.x,b.x),t),dec(add(mul(a.x,b.y),mul(a.y,b.x)),t));
	}
	friend cp operator*(cs cp &a,int b){return cp(mul(a.x,b),mul(a.y,b));}
	void operator+=(cs cp &b){Inc(x,b.x),Inc(y,b.y);}
	void operator-=(cs cp &b){Dec(x,b.x),Dec(y,b.y);}
	void operator*=(cs cp &b){*this=*this*b;}
	void operator*=(int b){Mul(x,b),Mul(y,b);}
};
inline cp inv(cs cp &a){
	int iv=inv(dec(add(mul(a.x,a.x),mul(a.y,a.y)),mul(a.x,a.y)));
	return cp(mul(dec(a.x,a.y),iv),mul(dec(0,a.y),iv));
}

cs int N=107,M=N*N;

int n,m;
cp a[N][N];

inline cp det(){
	cp ans(1);
	for(int re i=1;i<n;++i){
		int p;for(p=i;p<n;++p)if(a[p][i])break;
		if(p==n)return cp(0);
		if(p!=i){
			ans=cp(0)-ans;
			for(int re j=i;j<n;++j)std::swap(a[p][j],a[i][j]);
		}
		cp iv=inv(a[i][i]);ans*=a[i][i];
		for(int re j=i;j<n;++j)a[i][j]*=iv;
		assert(a[i][i].x==1);assert(a[i][i].y==0);
		for(int re j=i+1;j<n;++j)if(cp t=a[j][i])
		for(int re k=i;k<n;++k)a[j][k]-=t*a[i][k];
	}
	return ans;
}

cs cp omega[3]={cp(1,0),cp(0,1),cp(mod-1,mod-1)};
cs cp pw[3][3]={
	{omega[0],omega[0],omega[0]},
	{omega[0],omega[1],omega[2]},
	{omega[0],omega[2],omega[1]},
};

inline cp w1(cs cp &a){return cp(dec(0,a.y),dec(a.x,a.y));}
inline cp w2(cs cp &a){return cp(dec(a.y,a.x),dec(0,a.x));}

int u[M],v[M],w[M],c[M];

signed main(){
	freopen("sum.in","r",stdin);freopen("sum.out","w",stdout);
	n=gi(),m=gi();int mxc=0;
	for(int re i=1;i<=m;++i){
		u[i]=gi(),v[i]=gi();
		mxc=std::max(mxc,w[i]=gi());
	}
	int ans=0;
	for(int re p=1;p<=mxc;p*=3){
		for(int re i=1;i<=m;++i)c[i]=(w[i]/p)%3;
		static cp tp[3];
		for(int re i=0;i<3;++i){
			memset(a,0,sizeof a);
			for(int re j=1;j<=m;++j){
				int u=::u[j],v=::v[j];
				cp t=pw[i][c[j]];
				a[u][u]+=t;
				a[v][v]+=t;
				a[u][v]-=t;
				a[v][u]-=t;
			}
			tp[i]=det();
		}
		cp x=tp[0]+w2(tp[1])+w1(tp[2]);
		cp y=tp[0]+w1(tp[1])+w2(tp[2]);
		Inc(ans,mul(add(x.x,add(y.x,y.x)),p));
	}
	cout<<mul(ans,inv3)<<"\n";
	return 0;
}