6461. GDOI2020模拟02.05生成树（矩阵树及其扩展二维拉格朗日插值）

Posted 2021-03-07 gmh77

题目描述

给定一张 N 个点，M 条边的无向图，边有红、绿、蓝三种颜色，分别用 1,2,3 表示。

求这张图有多少生成树，满足绿色边数量不超过 x，蓝色边数量不超过 y，答案对10^9 + 7 取模。

1 ≤ N ≤ 40,1 ≤ M ≤ 10^5,1 ≤ ci ≤ 3

行列式

定义矩阵A的行列式det(A)或|A|

(|A|=sum_{排列p}{(-1)^{p的逆序对个数}prod{A_{i,p[i]}}})

行列式的性质

（(A)的转置矩阵(A^T)：把(A)的行列互换）

①(|A|=|A^T|)

②(|A||B|=|AB|)

见 https://www.zhihu.com/question/48497108

③(|A的某一行乘x|=x|A|)

④(|A的某行第j项=aj+bj|=|A的某行第j项=aj|+|A的某行第j项=bj|)

⑤交换相邻两行行列式变号

因为在每一种排列中符号都会改变

⑥交换任意两行行列式变号

两两依次交换的次数为2k+1，根据⑤可得必然变号

⑦有两行完全相同时行列式为0

否则交换两行后行列式改变

⑧把某一行乘k后加到另一行上行列式不变

相当于左乘了一个对角线为1，某个位置为k的矩阵， 这个矩阵的行列式显然是1

根据②可知不变

⑨上三角矩阵的行列式为对角线的积

其实可以不用化成上三角矩阵，只需要消元后n^2算逆序对即可

有了上面这些定理即可把矩阵高斯消元后求得行列式

（注意某行单独乘x时最终的结果要乘上x）

矩阵树定理

余子式：(M_{i,j})表示矩阵(A)去掉第i行第j列后剩余矩阵的行列式

基尔霍夫矩阵：最简单的形式即为 度数矩阵D-邻接矩阵C，其中度数矩阵(D_{i,i}=)i的度数

其实度数矩阵的真正形态是这样的：(D_{i,i}=sum_{i≠j}{C[i][j]})

所以当Ci,j不为1，甚至是一个多项式时也是成立的

矩阵树定理：图的生成树个数=基尔霍夫矩阵的任意余子式(M_{i,i})

我也不会证

二维拉格朗日插值

已知x=xi,y=yj时的点值(a_{x_iy_j})，构造(sum_{i,j}{a_{x_iy_j}f_{x_i}(x)}g_{x_i}(x))使得(x=x_i)时(f_{x_i}(x)=1)，否则=0，g同理

可以发现和一维的一模一样，(f_{x_i}(x)=prod_{j≠i}{frac{x-x_j}{x_i-x_j}})

题解

知道了上面的这题就是模板题了

设三种边边权为1,x,y，那么答案就是(x^Xy^Y)（XY为输入的数）的系数

枚举xy，求出行列式后插值即可

code

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdlib>
#include <cstdio>
#define fo(a,b,c) for (a=b; a<=c; a++)
#define fd(a,b,c) for (a=b; a>=c; a--)
#define mod 1000000007
#define Mod 1000000005
#define file
using namespace std;

long long A[41][41];
int a[41][41][3];
long long b[41][41];
long long c[41][41];
long long g[41];
long long ans[41][41];
int f[41];
int F[41];
int n,m,X,Y,i,j,k,l,x,y;
long long s,Ans;

long long qpower(long long a,int b)
{
    long long ans=1;
    
    while (b)
    {
        if (b&1)
        ans=ans*a%mod;
        
        a=a*a%mod;
        b>>=1;
    }
    
    return ans;
}

int main()
{
    freopen("tree.in","r",stdin);
    #ifdef file
    freopen("tree.out","w",stdout);
    #endif
    
    scanf("%d%d%d%d",&n,&m,&X,&Y);
    fo(i,1,m)
    scanf("%d%d%d",&j,&k,&l),--l,++a[j][k][l],++a[k][j][l];
    
    fo(x,0,n-1)
    {
        fo(y,0,n-1)
        {
            A[x][y]=1;
            
            fo(i,1,n)
            {
                b[i][i]=0;
                
                fo(j,1,n)
                if (i!=j)
                b[i][j]=-(a[i][j][0]+a[i][j][1]*x+a[i][j][2]*y),b[i][i]-=b[i][j];
            }
            
            memset(f,0,sizeof(f));
            fo(i,1,n-1)
            {
                fo(j,1,n-1)
                if (b[i][j])
                {
                    if (!f[j])
                    {
                        A[x][y]=A[x][y]*b[i][j]%mod;
                        
                        b[i][j]=qpower(b[i][j],Mod);
                        fo(k,j+1,n-1)
                        b[i][k]=b[i][k]*b[i][j]%mod;
                        b[i][j]=1;
                        
                        f[j]=i;
                        F[i]=j;
                        
                        fo(k,1,i-1)
                        if (F[k]>F[i])
                        A[x][y]=-A[x][y];
                        break;
                    }
                    else
                    {
                        fo(k,j+1,n-1)
                        b[i][k]=(b[i][k]-b[i][j]*b[f[j]][k])%mod;
                        b[i][j]=0;
                    }
                }
                
                if (j>n-1)
                {
                    A[x][y]=0;
                    break;
                }
            }
        }
    }
    
    fo(i,0,n-1)
    {
        c[i][0]=s=1;
        
        fo(j,0,n-1)
        if (i!=j)
        {
            memset(g,0,sizeof(g));
            s=s*(i-j)%mod;
            
            fo(k,0,n-2)
            {
                g[k]=(g[k]-c[i][k]*j)%mod;
                g[k+1]=(g[k+1]+c[i][k])%mod;
            }
            fo(k,0,n-1) c[i][k]=g[k];
        }
        
        s=qpower(s,Mod);
        fo(j,0,n-1)
        c[i][j]=c[i][j]*s%mod;
    }
    
    fo(i,0,n-1)
    {
        memset(g,0,sizeof(g));
        fo(j,0,n-1)
        {
            fo(k,0,n-1)
            g[k]=(g[k]+c[j][k]*A[i][j])%mod;
        }
        
        fo(j,0,n-1)
        {
            fo(k,0,n-1)
            ans[j][k]=(ans[j][k]+c[i][j]*g[k])%mod;
        }
    }
    
    fo(i,0,X)
    {
        fo(j,0,Y)
        Ans=(Ans+ans[i][j])%mod;
    }
    
    printf("%lld
",(Ans+mod)%mod);
    
    fclose(stdin);
    fclose(stdout);
    
    return 0;
}