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();}
posted @ 2020-07-24 20:14  RikukiIX  阅读(292)  评论(0编辑  收藏  举报