[BZOJ1005][HNOI2008]明明的烦恼

这一道题属于树的计数问题，我所知道的有三种方式处理这类问题，dp，基尔霍夫矩阵还有purfer编码。
这一题跟度数有关那么就用purfer编码啦。
purfer编码和树是一一对应的，如何把一棵树转化成编码。首先把树的度数为1的编号最小的点所连向的点放到编码的第一个，然后删除这个度为1的点，一直这样搞直到树的点数小于2为止。
一颗n的树的purfer编码只有n-2个数。
如何从purfer编码转化成一棵树，首先列出每个点的度数，即每个点在编码中出现的次数+1。选择编号最小的度数为1的点与编码的第一个点连边，把这两的点的度数都减1，然后再选编号最小的度数为1的点与编码的第二个点连边，依此类推。
本题无非就是要求满足条件的编码的个数。
每个点在编码中出现的次数为deg[i]-1，记这些次数的和为sum,有限制的点的个数为cnt。
那么有度数限制的点在purfer编码中本质不同的排列有：
$C^{sum}_{n-2-sum}*\frac{sum!}{\Pi_{i=1}^{cnt}(deg[i]-1)}$剩下n-2-sum个位置可以随便分给n-cnt个deg=-1的点，有$（n-cnt)^{(n-sum-2)}$种。
一共$C^{sum}_{n-2-sum}\times\frac{sum!}{\Pi_{i=1}^{cnt}(deg[i]-1)}\times（n-cnt)^{(n-sum-2)}$种。
可以化简，由于答案是整数，就不用写高精度除法了，直接分解质因数约分高精度乘就可以。

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<iostream>
#include<cassert>
#include<map>
#include<climits>
#include<vector>
#include<queue>
#include<cmath>
#define pb push_back
#define DB double
#define X first
#define Y second
#define MP make_pair
#define INF 0x3f3f3f3f
#define LL long long
#define pii pair<int,int>
#define pll pair<LL,LL>
#define DEBUG(...) fprintf(stderr,__VA_ARGS__)
using namespace std;
template<typename T>void Read(T& x)
{
    x=0;char c;int flag=0,sgn=1;
    while(c=getchar())
    {
        if(c=='-')sgn=-1;
        else if(c>='0'&&c<='9')x*=10,x+=c-'0',flag=1;
        else if(flag)break;
    }
    x*=sgn;
}
const int BASE=1e4,SIZE=1e4;
struct Bigint{
    int a[SIZE];
    void init(int num)
    {
        memset(a,0,sizeof(a));
        while(num)
            a[++a[0]]=num%BASE,num/=BASE;       
    }
    Bigint operator *(const Bigint& rhs)const{
        Bigint res;
        res.init(0);
        for(int i=1;i<=a[0];i++)
            for(int j=1;j<=rhs.a[0];j++)
            {
                res.a[i+j-1]+=a[i]*rhs.a[j];
                if(res.a[i+j-1]>BASE)res.a[i+j]+=res.a[i+j-1]/BASE,res.a[i+j-1]%=BASE;
            }
        for(int i=SIZE-1;i>=0;i--)
            if(res.a[i])
            {
                res.a[0]=i;
                break;
            }
        return res;
    }
    void prt()
    {
        printf("%d",a[a[0]]);
        for(int i=a[0]-1;i>=1;i--)
            printf("%04d",a[i]);
    }
};
int Cnt[5000],n,deg[5000],sum=0,pr[5000],flag[5000],tot=0;
void Euler()
{
    for(int i=2;i<5000;i++)
    {
        if(!flag[i])
            pr[tot++]=i,flag[i]=1;
        for(int j=0;j<tot&&pr[j]*i<5000;j++)
        {
            flag[i*pr[j]]=1;
            if(i%pr[j]==0)
                break;
        }
    }
}
void div(int k,vector<pii>& sol)
{
    for(int i=0;i<tot;i++)
    {
        int cnt=0;
        if(pr[tot]>sqrt(k))break;
        while(k%pr[i]==0)
            k/=pr[i],cnt++;
        if(cnt)
            sol.pb(MP(pr[i],cnt));
    }
    if(k!=1)sol.pb(MP(k,1));
}
Bigint ans,tmp,now;
int main()
{
    freopen("1005.in","r",stdin);
    freopen("1005.out","w",stdout);
    Read(n);
    Euler();
    int lover=0;
    for(int i=1;i<=n;i++)
    {
        Read(deg[i]);
        if(deg[i]!=-1)
            sum+=deg[i]-1,++lover;
    }
    if(sum>n-2)
    {
        puts("0");
        return 0;
    }
    vector<pii> s;
    for(int i=n-2;i>=n-2-sum+1;i--)
    {
        s.clear();
        div(i,s);
        for(int j=0;j<s.size();j++)
            Cnt[s[j].X]+=s[j].Y;
    }
    s.clear();
    div(n-lover,s);
    for(int j=0;j<s.size();j++)
            Cnt[s[j].X]+=s[j].Y*(n-2-sum);
    for(int i=1;i<=n;i++)
    {
        if(deg[i]==-1)continue;
        for(int k=deg[i]-1;k>=2;k--)
        {
            s.clear();
            div(k,s);
            for(int j=0;j<s.size();j++)
                Cnt[s[j].X]-=s[j].Y;
        }
    }

    ans.init(1);
    for(int i=1;i<5000;i++)
    {
        if(Cnt[i])
        {
            tmp.init(1);now.init(i);
            for(int j=1;j<=Cnt[i];j++)
                tmp=tmp*now;
            ans=ans*tmp;
        }
    }
    ans.prt();
}