树上M乘积路径数

此题是在某个的oj上看到的,题目标题和数据都不记得了,只能说个大概(标题瞎编的)

题目大意

给出\(n\)个点,可以任意选择生成树和\([1,m]\)的边权。由这些点连成一棵树,树中每一条边的权重介于\([1,m]\)给出两点\(a,b\)问有多少种树满足\(a,b\)之间的路径的权值乘积为\(m\)

数据范围

\(n<=10^{6}\)(应该是),\(m<=10^{18}\)

题目解析

这里有道简化版的类似题目:CodeForces #539 (Div.2) F:Sasha and Interesting Fact from Graph Theory

首先此题需要用到凯利公式(原题提示中已给出):$$设T_{n,k}为n个有标号点分成k个连通块的的分割种数,其中节点1,2,3,…k属于不同的联通块,则有
T_{n,k}=k\ n^{n-k-1}$$

与上面简化版类似,假设给定的\(a,b\)之间假设存在条路径,总共有\(i\)条边,\(i+1\)个点,那么除了\(a,b\)选择剩下的节点方案数为\(C_{n-2}^{i-1}\),乘以对点进行全排列\(A_{i-1}\)

满足ab路径的边权条件后,其他边权的分配方案为\(m^{n-1-i}\)

然后将剩下的\(n-i-1\)个节点接入\(a,b\)的路径,除去\(a,b\)路径上的边,就转化为了将\(n\)个点划分到\(i+1\)个联通块(且\(a,b\)路径上的i+1个点属于不同联通块),这就可以使用凯利公式\(T_{n,i+1}=(i+1)*n^{n-i-2}\)

接下来就是计算\(i+1\)条边权乘积为m的方案数。如果是简化版题目中的和为\(m\),则只要使用隔板法就可以用一个组合数表示,然而乘积的情况下,就等对原题目进行质因数分解,转化为质因数的指数幂形式,在对每个质因数的指数分别进行隔板法再结果相乘(而且每个边代表的指数是可以为0的)。假设质因数p的指数幂为k,则隔板法方案数为\(C_{k+i-1}^{i-1}\)

但是本题m非常大,常规质因数分解算法\(O(\sqrt{m})\)的质因数分解算法无法使用,这里想到的是使用大数因子分解算法Pollard Rho结合Miller Rabin提取素数因子。

如果提示中不给出凯利公式,那应该就是一道高难度的题了(给出了难度也不低就是)

代码

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<ctime>
#include<map>
using namespace std;
typedef long long ll;
const int MAXN = 1000005 ;
const int MOD=1e9+7;

///logn快速乘法
ll mulEx(ll a , ll b , ll Mod) {
    if(!a) return 0 ;
    ll ans(0) ;
    while(b)
    {
		if(b & 1) ans = (ans + a) % Mod;
		a <<= 1 ;
		a %= Mod ;
		b >>= 1 ;
    }
    return ans ;
}
 
///快速幂
ll powEx(ll base , ll n , ll Mod)
{
    ll ans(1) ;
    while(n)
    {
        if(n & 1) ans = mulEx(ans , base , Mod) ;
        base = mulEx(base , base , Mod) ;
        n >>= 1 ;
    }
    return ans ;
}
 
bool check(ll a , ll d , ll n)
{
    if(n == a) return true ;
    while(~d & 1) d >>= 1 ;
    ll t = powEx(a , d , n) ;
    while(d < n - 1 && t != 1 && t != n - 1)
    {
        t = mulEx(t , t , n) ;
        d <<= 1 ;
    }
    return (d & 1) || t == n - 1 ;
}

///判断大数是否是质数,miller_rabin
bool isP(ll n)
{ 
    if(n == 2) return true ;
    if(n < 2 || 0 == (n & 1)) return false ;
    static int p[5] = {2 , 3 , 7 , 61 , 24251} ;
    for(int i = 0 ; i < 5 ; ++ i) if(!check(p[i] , n - 1 , n)) return false ;
    return true ;
}
 
ll gcd(ll a , ll b)
{
    if(a < 0) return gcd(-a , b) ;
    return b ? gcd(b , a - b * (a / b)) : a ;
}
 
///大数分解质因数
ll Pollard_rho(ll n , ll c)
{
    ll i = 1 , k = 2 , x = rand() % n , y = x ;
    while(true)
    {
        x = (mulEx(x , x , n) + c) % n ;
        ll d = gcd(y - x , n) ;
        if(d != 1 && d != n) return d ;
        if(y == x) return n ;
        if(++ i == k)
        {
            y = x ;
            k <<= 1 ;
        }
    }
}
 
ll Fac[MAXN] , factCnt ;
///Fac存的是质因子,大小不一定按照顺序,有重复
map<ll , ll> factMap ;
///遍历map的first表示因子,second表示次数

void factorization(ll n)
{
    if(isP(n))
    {
        if (factMap.count(n)==0){
            Fac[factCnt++] = n ;
            factMap[n]=0;
        }
        return ;
    }
    ll p(n) ;
    while(p >= n) p = Pollard_rho(p , rand() % (n - 1) + 1) ;
    factorization(p) ;
    factorization(n / p) ;
}
 
 
void getFactor(ll x)
{
    if (isP(x)){
        Fac[factCnt++]=x;
        factMap[x]=0;
        return;
    }
	srand(time(0)) ;
	factCnt = 0 ;
	factMap.clear() ;
	factorization(x) ;
}
ll fac[2000005],ifac[2000005],nexp[1000005],mexp[1000005];
void init(ll n,ll m){
    ifac[1]=1;ifac[0]=1;
    fac[0]=1;fac[1]=1;
    for (ll i=2;i<=2*n;i++){
        fac[i]=fac[i-1]*i%MOD;
        ifac[i]=1ll*(MOD-MOD/i)*ifac[MOD%i]%MOD;
    }
    for (ll i=1;i<=2*n;i++){
        ifac[i]=ifac[i-1]*ifac[i]%MOD;
    }
    nexp[0]=1;nexp[1]=n;
    mexp[0]=1;mexp[1]=m%MOD;
    for (int i=2;i<=n;i++){
        nexp[i]=nexp[i-1]*n%MOD;
        mexp[i]=mexp[i-1]*(m%MOD)%MOD;
    }
}
//x>y
ll getC(int x,int y){
    return fac[x]*ifac[y]%MOD*ifac[x-y]%MOD;
}
int main(){
    int a,b;
    ll m,n;
    cin>>n>>m>>a>>b;
    init(n,m);
    ll invn=powEx(n,MOD-2,MOD);

    if (m==1){
        ll ans=0;
        for (int i=1;i<=n-1;i++){
            ll sums=1;
            sums=sums*fac[i-1]%MOD;
            sums=sums*getC(n-2,i-1)%MOD;
            if (n-i-2>=0)
                sums=sums*(i+1)%MOD*nexp[n-i-2]%MOD;
            else sums=sums*(i+1)%MOD*powEx(n,MOD-2,MOD)%MOD;
            ans=(ans+sums)%MOD;
        }
        cout<<ans;
        return 0;
    }
    getFactor(m);
    ll v=m;
    for (int i=0;i<factCnt;i++){
        while (v%Fac[i]==0){
            factMap[Fac[i]]++;
            v/=Fac[i];
        }
    }
    ll ans=0;
    for (int i=1;i<=n-1;i++){
        ll sums=1;
        for (int j=0;j<factCnt;j++){
            int tj=factMap[Fac[j]];
            sums=sums*getC(tj+i-1,i-1)%MOD;
        }
        sums=sums*fac[i-1]%MOD;
        sums=sums*mexp[n-1-i]%MOD;
        sums=sums*getC(n-2,i-1)%MOD;
        if (n-i-2>=0)
            sums=sums*(i+1)%MOD*nexp[n-i-2]%MOD;
        else sums=sums*(i+1)%MOD*invn;
        ans=(ans+sums)%MOD;
    }
    cout<<ans;
}```
posted @ 2022-07-26 15:48  Z-Y-Y-S  阅读(77)  评论(0编辑  收藏  举报