1. 程式人生 > 實用技巧 >loj6271 「長樂集訓 2017 Day10」生成樹求和 加強版(矩陣樹定理,迴圈卷積)

loj6271 「長樂集訓 2017 Day10」生成樹求和 加強版(矩陣樹定理,迴圈卷積)

loj6271 「長樂集訓 2017 Day10」生成樹求和 加強版(矩陣樹定理,迴圈卷積)

loj

題解時間

首先想到先分開三進位制下每一位,然後每一位分別求結果為0,1,2的樹的個數。

然後考慮矩陣樹求出來的是樹的邊權之積的和,而我們要求樹的邊權的不進位三進位制和的和。

由於矩陣樹求出來的是樹的邊權之積的和,考慮暴力上生成函式求解迴圈卷積,結果就是 $ c $ 的項的係數和。

但很明顯生成函式暴力算是沒得整的。

所以我們想到了利用單位根實現的k進位制FWT。

很幸運的 $ \omega_{ 3 }^{ 2 } = - \omega_{ 3 } -1 $ 。

這樣一來數值可以直接用一個 $x + y \omega_{ 3 } $ 來表示。

求出每一位的點值之後直接IFWT回去。

但是我不會(悲)

程式碼

為什麼會有檔案讀寫啊,正巧在loj剛咕完恢復過來的時候交,T了好多次以為loj評測姬也咕了呢。。。

不就是你自己不好好讀題。

#include<bits/stdc++.h>
using namespace std;
typedef long long lint;
struct pat{int x,y;pat(int x=0,int y=0):x(x),y(y){}bool operator<(const pat &p)const{return x==p.x?y<p.y:x<p.x;}};
template<typename TP>inline void read(TP &tar)
{
	TP ret=0,f=1;char ch=getchar();
	while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
	while(ch>='0'&&ch<='9'){ret=ret*10+(ch-'0');ch=getchar();}
	tar=ret*f;
}
template<typename TP,typename... Args>inline void read(TP& t,Args&... args){read(t),read(args...);}
namespace RKK
{
const int N=110,M=5011,C=19683;
const int mo=1000000007,inv3=333333336;
void doadd(int &a,int b){if((a+=b)>=mo) a-=mo;}int add(int a,int b){return (a+=b)>=mo?a-mo:a;}
void dodec(int &a,int b){doadd(a,mo-b);}int dec(int a,int b){return add(a,mo-b);}
void domul(int &a,int b){a=1ll*a*b%mo;}int mul(int a,int b){return 1ll*a*b%mo;}
int fpow(int a,int p){int ret=1;while(p){if(p&1) domul(ret,a);domul(a,a),p>>=1;}return ret;}
struct cp
{
	int x,y;//x+yw
	cp(const int &x=0,const int &y=0):x(x),y(y){}
	operator bool(){return x||y;}
	cp inv()const{int i=fpow(dec(add(mul(x,x),mul(y,y)),mul(x,y)),mo-2);return cp(mul(dec(x,y),i),mul(dec(0,y),i));}
	cp operator+(const cp &p)const{return cp(add(x,p.x),add(y,p.y));}
	cp operator-(const cp &p)const{return cp(dec(x,p.x),dec(y,p.y));}
	cp operator*(const cp &p)const{return cp(dec(mul(x,p.x),mul(y,p.y)),dec(add(mul(x,p.y),mul(p.x,y)),mul(y,p.y)));}
	void operator+=(const cp &p){*this=*this+p;}
	void operator-=(const cp &p){*this=*this-p;}
	void operator*=(const cp &p){*this=*this*p;}
};
const cp om[]={cp(1,0),cp(0,1),cp(mo-1,mo-1)};
int n,m,ans,ex[M],ey[M],ew[M];
cp a[N][N];
cp calc()
{
	cp ret=cp(1,0);
	for(int l=1;l<n;l++)
	{
		int e=0;for(int i=l;i<n;i++)if(a[i][l]){e=i;break;}
		if(!e) return cp(0,0);
		if(e!=l){ret=cp(0,0)-ret;for(int j=1;j<n;j++) swap(a[l][j],a[e][j]);}
		cp inv=a[l][l].inv();
		for(int i=l+1;i<n;i++)
		{
			cp k=a[i][l]*inv;
			for(int j=l;j<n;j++) a[i][j]-=k*a[l][j];
		}
		ret*=a[l][l];
	}return ret;
}
int main()
{
	freopen("sum.in","r",stdin),freopen("sum.out","w",stdout);
	read(n,m);for(int i=1;i<=m;i++) read(ex[i],ey[i],ew[i]);
	for(int p=1;p<C;p*=3)
	{
		static cp at[3];
		static cp w,x,y;
		for(int i=0;i<3;i++)
		{
			memset(a,0,sizeof(a));
			for(int j=1;j<=m;j++)
			{
				w=om[(i*((ew[j]/p)%3))%3];
				a[ex[j]][ex[j]]+=w,a[ey[j]][ey[j]]+=w;
				a[ex[j]][ey[j]]-=w,a[ey[j]][ex[j]]-=w;
			}at[i]=calc();
		}
		x=at[0]+at[1]*om[2]+at[2]*om[1],y=at[0]+at[1]*om[1]+at[2]*om[2];
		doadd(ans,mul(p,add(x.x,add(y.x,y.x))));
	}
	domul(ans,inv3);
	printf("%d\n",ans);
	return 0;
}
}
int main(){return RKK::main();}