4.9 省选模拟赛 生成树求和 变元矩阵树定理 生成函数 iDFT 插值法
有同学在loj上找到了加强版 所以这道题是可以交的。LINK:生成树求和 加强版
对于30分 爆搜 可实际上我爆搜只过了25分 有同学使用按秩合并并茶几的及时剪枝通过了30分。
const int MAXN=45;
int n,m;
struct wy
{
int x,y,z;
}t[MAXN];
int w[MAXN],f[MAXN];ll ans;
inline int getfather(int x){return x==f[x]?x:getfather(f[x]);}
inline int add(int x,int y)
{
int cnt=0;
int p=1;
while(x||y)
{
int w1=x%3;
int w2=y%3;
int cc=(w1+w2)%3;
cnt+=cc*p;
p=p*3;x/=3;y/=3;
}
return cnt;
}
inline void dfs(int x,int v)
{
if(m-x+1<n-v)return;
if(v==n)
{
int cnt=0;
rep(1,n,i)f[i]=i;
rep(1,n-1,i)
{
int xx=getfather(t[w[i]].x);
int yy=getfather(t[w[i]].y);
if(xx==yy)return;
f[xx]=yy;
cnt=add(cnt,t[w[i]].z);
}
ans=(ans+cnt)%mod;
return;
}
w[v]=x;
dfs(x+1,v+1);
dfs(x+1,v);
}
int main()
{
freopen("sum.in","r",stdin);
freopen("sum.out","w",stdout);
get(n);get(m);
rep(1,m,i)
{
int get(x);int get(y);int get(z);
t[i]=(wy){x,y,z};
}
if(m<=35){dfs(1,1);putl(ans);}
return 0;
}
之所以这样做是因为虽然这道题很容易想到矩阵树定理,按位处理 但是依然无法做出来。
矩阵树除了能求方案还可以求边上的权值积的和,即变元矩阵树定理 具体理解可以把边权想象成有边权条重边。
但是对于这道题来说还是不行 求的是生成树边权和的和 把积变和只有是在指数上了。
而且还有三进制的加法 非常的不可做 但是可以按位处理。
把边权变到指数上只能是 生成函数了 列出生成函数 每一条边的权值用生成函数表示。
利用高斯消元解行列式 最后求出的多项式就好了。
而直接使用多项式无法方便的求出行列式 中间由于矩阵树定理还无法取模 所以有多项式乘法除法啥的 所以变的更不可做了。
但是考虑最后的次数 其实并不大最多2n次 可以取2n+1个数带入 然后利用插值求解原多项式。
复杂度\(n^4log\)
为了通过加强版 我们考虑一种比较简单的方式 如果能找到\(x^3=1(mod P)\)
那么我么最后的多项式就只有三项且是正确的了 考虑到复平面单位根。
由于最后有三项 所以 把复平面平分三份的复平面单位根即可。
考虑最后的插值 利用复平面单位根的相加相乘操作可以轻松的插值出来。不需要再次高斯消元或者拉格朗日插值。
非常的巧妙。
const ll MAXN=100010,inv2=(mod+1)>>1,inv3=(mod+1)/3,g3=82062379,maxn=110;
inline ll ksm(ll b,ll p,ll cc){ll cnt=1;while(p){if(p&1)cnt=cnt*b%cc;p=p>>1;b=b*b%cc;}return cnt;}
struct wy
{
ll a,b;
inline wy(ll aa=0,ll bb=0){a=aa;b=bb;}
inline wy operator +(wy w){return wy((a+w.a)%mod,(b+w.b)%mod);}
inline wy operator -(wy w){return wy((a-w.a+mod)%mod,(b-w.b+mod)%mod);}
inline wy operator *(ll y){return wy(a*y%mod,b*y%mod);}
inline wy operator *(wy y){return wy((a*y.a%mod-b*y.b%mod+mod)%mod,(a*y.b+b*y.a)%mod);}
inline wy operator /(ll y)
{
ll inv=ksm(y,mod-2,mod);
return wy(a*inv%mod,b*inv%mod);
}
inline wy operator /(wy y)
{
return wy(a,b)*wy(y.a,mod-y.b)/((y.a*y.a+y.b*y.b)%mod);
}
}a[maxn][maxn],cs[3],inv,w[3],invw[3];
int pow3[maxn],eu[maxn*maxn],ev[maxn*maxn],ec[maxn*maxn];
ll n,m;
inline wy det()
{
wy ret=wy(1,0);
rep(1,n-1,i)
{
ll p=i;
rep(i,n-1,j)if(a[j][i].a||a[j][i].b){p=j;break;}
if(i!=p)
{
rep(1,n-1,k)swap(a[i][k],a[p][k]);
ret=ret*(mod-1);
}
rep(i+1,n-1,j)
{
inv=a[j][i]/a[i][i];
rep(i,n-1,k)a[j][k]=a[j][k]-inv*a[i][k];
}
ret=ret*a[i][i];
}
return ret;
}
inline void idft(wy *p)
{
wy c[3];
c[0]=p[0];c[1]=p[1];c[2]=p[2];
p[0]=c[0]+c[1]+c[2];
p[1]=c[0]+invw[1]*c[1]+c[2]*invw[2];
p[2]=c[0]+invw[2]*c[1]+c[2]*invw[1];
rep(0,2,i)p[i]=p[i]*inv3;
}
inline ll calc(ll p)
{
wy w0,w1,w2,c;
rep(0,2,i)
{
memset(a,0,sizeof(a));
w0=wy(1,0);
w1=w[i];
w2=w[(i+i)%3];
rep(1,m,j)
{
ll k=ec[j]/pow3[p]%3;
c=(k==1?w1:k==2?w2:w0);
a[eu[j]][eu[j]]=a[eu[j]][eu[j]]+c;
a[ev[j]][ev[j]]=a[ev[j]][ev[j]]+c;
a[eu[j]][ev[j]]=a[eu[j]][ev[j]]-c;
a[ev[j]][eu[j]]=a[ev[j]][eu[j]]-c;
}
cs[i]=det();
}
idft(cs);
return (cs[1].a+cs[2].a*2)*pow3[p]%mod;
}
signed main()
{
freopen("1.in","r",stdin);
//freopen("sum.out","w",stdout);
get(n);get(m);ll ans=0;
w[0]=wy(1,0);w[1]=wy(mod-inv2,g3*inv2%mod);
w[2]=wy(mod-inv2,mod-g3*inv2%mod);
invw[0]=w[0]/w[0];
invw[1]=w[0]/w[1];
invw[2]=w[0]/w[2];
rep(1,m,i)get(eu[i]),get(ev[i]),get(ec[i]);
pow3[0]=1;
rep(1,10,i)pow3[i]=pow3[i-1]*3;
rep(0,10,i)ans=(ans+calc(i))%mod;
putl((ans+mod)%mod);return 0;
}