BZOJ4271：chemistry 化学

【题目描述】

烷烃，即饱和烃（saturated group)，是只有碳碳单键和碳氢键的链烃，是最简单的一类有机化合物。

化学上，同分异构体是一种有相同化学式，有同样的化学键而有不同的原子排列的化合物。简单地说，化合物具有相同分子式，但具有不同结构的现象，叫做同分异构现象；具有相同分子式而结构不同的化合物互为同分异构体。很多同分异构体有相似的性质。

本题要求n烷的同分异构体个数。例如，丁烷（四烷？）有两个：正丁烷，异丁烷。

【输入】
第一行：n，含义见题意。

【输出】
第一行：答案。

【输入样例】
4
【输出样例】
2
【提示】
对于100%的数据，1≤n≤500。

##----------前言----------
####某天大课间：
watson：你们学没学烷烃的同分异构体啊，我发现烷烃的同分异构体就是每一个点度数不超过4的无根树个数嘛，貌似可以用OI姿势DP？？

我：那就是n个点度数不超过m的无根树个数？喵喵喵？我不会啊，这玩意怎么计数，怎么DP？……

####数学课上：
诶，对于n个点度数不超过m的有根树，貌似有一种O(n^(m - 1))的算法，考虑dp[i]表示有i个节点根节点度数不超过3的方案数，那么dp[i]可以从各种各样的子树中任选3个/2个/1个接在上面，这个大力枚举一下就可以有一个指数级算法啦~

然后呢？然后没然后了QAQ

####机房中：
嗯？__debug大佬貌似写过一篇关于烷烃的[博客](http://debug18.com/posts/calculate-the-number-of-structural-isomers-for-alkanes/)？前去学习之后仍然一脸懵逼……

####第二天语文课：
琢磨了一下__debug大佬的转移方程，突然发现有道理！用$f(i,j,k)$表示有$i$个点，最大子树为$j$，根节点度数为$k$的有根树有多少种。
dp方程就是：
$$f(i, j, k) += \sum^{n}_{s=1}\sum^{\frac{i}{j}}_{x=1}f(i-j\times x, s, k - x)\cdot \binom {\sum^{j}_{a=1}\sum^{m-1}_{b=0}f(j,a,b)+x-1} x $$

思路：考虑新添一个根节点，枚举按大小顺序增添子树，这样容易限定方案是唯一的。那么我们就可以枚举最大子树的大小$j$，再枚举大小为$j$的子树个数$x$，这样我们相当于要从以前所有大小为$j$的$\sum^{j}_{a=1}\sum^{m-1}_{b=0}f(j,a,b)$种子树(每种都有无限个，并且度数严格小于$m$，否则连接后度数会大于$m$)中任选出$x$个来接在当前的根节点上。就可以用组合数计数。

容易知道下面这个部分可以提前算，相当于降了一阶：
$$\sum^{j}_{a=1}\sum^{m-1}_{b=0}f(j,a,b)$$

考虑优化这个dp，发现可以省去枚举$s$，我们只需要对于$f(i,j,k)$的每一个$(i, k)$下标处理出在$j$这一维上的前缀和即可。

但是实际上我们还可以类似背包的思想，发现如果规定$s$的顺序转移，后面的只会从前面的转移，可以类似背包优化地规定最外层枚举$s$，这样就可以节省空间。

那么无根树呢？其实很简单，如果规定有根树的根为重心，那么无根树和有根树一一对应！所以只需要规定最大子树大小不超过总点数的一半即可！(但是有两个重心需要特殊处理，我们就可以先处理子树严格小于总大小的一半的，最后再枚举合并两棵相同的子树，加上这一部分答案)

结果高精度被卡常了……

```
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#define ll long long
#define For(i, a, b) for(register int i = a; i <= b; ++i)
using namespace std;
const int maxn = 5e2 + 5;
namespace BIGINT {
    const int POW = 10000, ML = 530, st = 4;
    char s[ML]; int len, tmp, tmp2;
    const int TP[] = {1, 10, 100, 1000, 10000, 100000, 1000000, 10000000, 100000000};
    #define For(i, a, b) for(register int i = a; i <= b; ++i)
    #define Down(i, a, b) for(register int i = a; i >= b; --i)
    struct INT1024 {
        int a[ML], l;
        INT1024() {}
        void clear() {memset(a, 0, sizeof(int) * (l + 1)), l = 0;}
        INT1024(int x) {
            clear();
            while(x) a[l++] = x % POW, x /= POW;
        }
        void init(int x) {
            clear();
            while(x) a[l++] = x % POW, x /= POW;
        }
        void put() {
            printf("%d", a[l - 1]);
            Down(i, l - 2, 0) printf("%04d", a[i]);
        }
        void read() {
            scanf("%s", s), len = strlen(s);
            l = tmp2 = tmp = 0;
            Down(i, len - 1, 0) {
                tmp = tmp + (s[i] - '0') * TP[tmp2++];
                if(tmp2 == st) a[l++] = tmp, tmp2 = tmp = 0;
            }
            if(tmp) a[l++] = tmp;
        }
    } ;
    int mid[ML << 1];
    INT1024 operator + (INT1024 a, const INT1024 &b) {
        a.l = max(a.l, b.l);
        For(i, 0, a.l - 1) a.a[i] += b.a[i];
        For(i, 0, a.l - 1) {
            a.a[i + 1] += a.a[i] / POW;
            a.a[i] %= POW;
        }
        while(a.a[a.l]) 
            a.a[a.l + 1] += a.a[a.l] / POW, a.a[(a.l++)] %= POW;
        return a;
    }
    INT1024 operator + (INT1024 a, const int &b) {
        a.a[0] += b;
        For(i, 0, a.l - 1) {
            a.a[i + 1] += a.a[i] / POW;
            a.a[i] %= POW;
        }
        while(a.a[a.l]) 
            a.a[a.l + 1] += a.a[a.l] / POW, a.a[(a.l++)] %= POW;
        return a;
    }
    void operator += (INT1024 & a, const int &b) {
        a.a[0] += b;
        For(i, 0, a.l - 1) {
            if(a.a[i] < POW) return ;
            a.a[i + 1] += a.a[i] / POW;
            a.a[i] %= POW;
        }
        while(a.a[a.l]) 
            a.a[a.l + 1] += a.a[a.l] / POW, a.a[(a.l++)] %= POW;
    }
    INT1024 operator - (INT1024 a, const INT1024 &b) {
        a.l = max(a.l, b.l);
        For(i, 0, a.l - 1) a.a[i] -= b.a[i];
        For(i, 0, a.l - 2) {
            if(a.a[i] < 0) a.a[i + 1] -= 1, a.a[i] += POW;
        }
        while(!a.a[a.l - 1]) --a.l;
        return a;
    }
    INT1024 operator - (INT1024 a, const int &b) {
        a.a[0] -= b;
        For(i, 0, a.l - 2) {
            if(a.a[i] < 0) a.a[i + 1] -= 1, a.a[i] += POW;
            else return a;
        }
        while(!a.a[a.l - 1]) --a.l;
        return a;
    }
    void operator -= (INT1024 & a, const int &b) {
        a.a[0] -= b;
        For(i, 0, a.l - 2) {
            if(a.a[i] < 0) a.a[i + 1] -= 1, a.a[i] += POW;
            else return ;
        }
        while(!a.a[a.l - 1]) --a.l;
    }
    long long tml;
    INT1024 operator / (INT1024 a, const int &b) {
        tml = 0;
        Down(i, a.l - 1, 0) {
            tml = tml * POW + a.a[i];
            a.a[i] = tml / b, tml %= b;
        }
        while(!a.a[a.l - 1]) --a.l;
        return a;
    }
    void operator /= (INT1024 & a, const int &b) {
        tml = 0;
        Down(i, a.l - 1, 0) {
            tml = tml * POW + a.a[i];
            a.a[i] = tml / b, tml %= b;
        }
        while(!a.a[a.l - 1]) --a.l;
    }
    INT1024 operator * (INT1024 a, const INT1024 &b) {
        int len = a.l + b.l;
        memset(mid, 0, sizeof(int) * (len + 1));
        For(i, 0, a.l - 1) For(j, 0, b.l - 1)
            mid[i + j] += a.a[i] * b.a[j];
        For(i, 0, len - 1) mid[i + 1] += mid[i] / POW, mid[i] %= POW;
        while(!mid[len - 1]) --len;
        a.l = len;
        For(i, 0, len - 1) a.a[i] = mid[i];
        return a;
    }
}
using namespace BIGINT;
int n, m = 4;
INT1024 dp[maxn][5], sum, ans;
INT1024 t;
void C(INT1024 a, const int & b) {
    t.init(1);
    For(i, 1, b) {
        t = t * a;
        a -= 1; 
    }
    For(i, 2, b) 
        t /= i;
}
 
int main() {
    dp[1][0].a[0] = 1, dp[1][0].l = 1;
    scanf("%d", &n);
    For(s, 1, (n - 1) >> 1) {
        sum.clear();
        For(i, 0, m - 1) sum = sum + dp[s][i];
        for(register int i = n; i >= s; --i)
        For(j, 1, m) For(k, 1, j){
            if(i <= s * k) break;
            C(sum + k - 1, k);
            dp[i][j] = dp[i][j] + dp[i - s * k][j - k] * t;
        }
        /*for(register int i = 1; i <= n; ++i) {
            cerr << i << " : " << endl;
            for(register int j = 0; j <= m; ++j) {
                cerr << j << " : ";dp[i][j].put();cerr << "   ";
            }
            cerr << endl;
        }
        cerr << "++++++++++++++++++++++++" << endl; */
    }/*
    for(register int i = 1; i <= n; ++i) {
        cerr << i << " : " << endl;
        for(register int j = 0; j <= m; ++j) {
            cerr << j << " : " << dp[i][j] << "   ";
        }
        cerr << endl;
    }*/
    sum.clear(), ans.clear();
    For(i, 0, m) ans = ans + dp[n][i];
    if(n % 2 == 0) {
        For(i, 0, m - 1) sum = sum + dp[n / 2][i];
        C(sum + 1, 2);
        ans = ans + t;
    }
    ans.put();
    return 0;
}
```

花絮：2018烷的同分异构体个数：
```
55646884355336227932095012099708927190751318887053748208624788192905565244687147018097703534576491556446670726618274518583887744093476349899078747009604647181967553569481656870822295818863816000570944280469165114459352037189555619717330057816790655770097784063961415669760001686183104700695154762061638836314043559145668218631221097931994202798660096093169294169668331552533138269188512967579659005400328023092580578134113122702506827421701352281675214520954345089899733919186624919968221468961203233762350330410282331847050554015302064618344144243699950369940870475658984995543080988345220996052003578647252855031917076526359272769022104927991053639291275634147359331232576638846190109360429323426136901509008532676726419991378059265942576470372218260653372116989769754186483899691205954666081060703680385215644123559299044261871071054408913494675483051091618406162578274866818812370795906590092140

```

Rockdu's Blog

导航

最近发表

大佬们的友链