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();}

巴特西

loj6271 「长乐集训 2017 Day10」生成树求和加强版(矩阵树定理，循环卷积)

loj6271 「长乐集训 2017 Day10」生成树求和加强版(矩阵树定理，循环卷积)

题解时间

代码

最新文章

热门文章

巴特西

loj6271 「长乐集训 2017 Day10」生成树求和 加强版(矩阵树定理，循环卷积)

loj6271 「长乐集训 2017 Day10」生成树求和 加强版(矩阵树定理，循环卷积)

题解时间

代码

最新文章

热门文章

loj6271 「长乐集训 2017 Day10」生成树求和加强版(矩阵树定理，循环卷积)

loj6271 「长乐集训 2017 Day10」生成树求和加强版(矩阵树定理，循环卷积)