[POJ-1322]Chocolate(生成函数)

题目

口袋里有$c$种颜色的巧克力。从中随便取出一个。若桌上有与它颜色相同的巧克力，则两个都拿走；否则将它放在桌子上。求拿出$n$个后桌子上剩$m$个的概率。
数据范围：$1\le c\le 100,n,m\le 1000000$

思路

1. 易得：n,m奇偶性相同
2. $m\le c,n$

设计

设有$a$种巧克力取了奇数次，则有$c-a$种取了偶数次。
考虑到，若是使用普通生成函数的话，因为结果有除掉$k!$，所以是不好做的，这个时候，我们就选择使用指数型生成函数(EGF)。
奇数次的生成函数：
$$J(x)=\frac{x}{1!}+\frac{x^3}{3!}+\frac{x^5}{5!}+...$$
偶数次的生成函数：
$$O(x)=1+\frac{x^2}{2!}+\frac{x^4}{4!}+...$$
设答案多项式为$H(x)$，有：
$$H(x)=C_c^a\times J(x{)}^a\times O(x{)}^{c-a}$$
$H(x)$的第$n$项前的系数就是合法的组合的方案数

转化

设$f(x)=e^x$，对$f(x)$在零点处进行泰勒展开。
$$f(x)=f(x_0)+\frac{f^{\prime }(0)}{1!}x+\frac{f^{\prime \prime }(0)}{2!}{x}^2+...$$
因为$e^x=e^{\prime x}$
所以$f(x)=f^{\prime }(x)$
上式可化为
$$f(x)=1+\frac{x}{1!}+\frac{x^2}{2!}+\frac{x^3}{3!}+\frac{x^4}{4!}+...$$
$$f(-x)=1-\frac{x}{1!}+\frac{x^2}{2!}-\frac{x^3}{3!}+\frac{x^4}{4!}+...$$
所以
$$e^x=1+\frac{x}{1!}+\frac{x^2}{2!}+\frac{x^3}{3!}+\frac{x^4}{4!}+...$$
$$e^{-x}=1-\frac{x}{1!}+\frac{x^2}{2!}-\frac{x^3}{3!}+\frac{x^4}{4!}+...$$
$$e^x-e^{-x}=2\frac{x}{1!}+2\frac{x^3}{3!}+...2\frac{x^5}{5!}$$
$$e^x+e^{-x}=2+2\frac{x^2}{2!}+2\frac{x^4}{4!}+...$$
所以奇数次的生成函数的闭形式为：
$$\frac{e^x-e^{-x}}{2}$$
偶数次的生成函数的闭形式为：
$$\frac{e^x+e^{-x}}{2}$$

化简

由上可得:
$$H(x)=C_c^a\times J(x{)}^a\times O(x{)}^{c-a}=(\frac{e^x-e^{-x}}{2}{)}^a\times (\frac{e^x+e^{-x}}{2}{)}^{c-a}$$
进行二项式展开：
$$H(x)=2^{-c}\times C_c^a\times \sum _{i=0}^a{-1}^i\times C_a^i\times e^{a-2i}x\times \sum _{j=0}^{c-a}{C}_{c-m}^j\times e^{c-a-2j}$$
写成卷积：
$$H(x)=2^{-c}\times C_c^a\times \sum _{i=0}^a\sum _{j=0}^{c-a}{-1}^i\times C_a^i\times C_{c-m}^j\times e^{c-2i-2j}x$$
把$e$换掉
$$H(x)=2^{-c}\times C_c^a\times \sum _{i=0}^a\sum _{j=0}^{c-a}{-1}^i\times C_a^i\times C_{c-m}^j\times \sum _{k=0}^{+\infty }\frac{((c-2i-2j)x{)}^k}{k!}$$
所以
$$h_n=2^{-c}\times C_c^a\times \sum _{i=0}^a\sum _{j=0}^{c-a}{-1}^i\times C_a^i\times C_{c-m}^j\times \frac{(c-2i-2j{)}^n}{n!}$$

得到答案

上面的$h_n$是合法的组合的方案数，因为颜色是有顺序的，所以还要乘上一个$n!$。而总方案数为$c^n$
所以最后的$Ans$为：
$$\frac{h_n\ast n!}{c^n}$$
因为$1\le c\le 100$,$2^c$会乘爆，只能用$double$乘一下，用精度换数值正确。
把除放在最后，避免精度误差。

Code:

#include <cmath>
#include <queue>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#define DB double
#define MAXN 105
#define LL long long
#define Int register int
using namespace std;
inline int read()
{
	int x = 0, f = 1;
	char s = getchar();
	while (s < '0' || s > '9')
	{
		if (s == '-')
			f = -1;
		s = getchar();
	}
	while (s >= '0' && s <= '9')
	{
		x = (x << 3) + (x << 1) + (s ^ 48);
		s = getchar();
	}
	return x * f;
}
inline int Max(int x,int y)
{
	return x > y ? x : y;
}
inline int Min(int x,int y)
{
	return x < y ? x : y;
}
inline DB QuickPow(DB x,int y)
{
	DB Res = 1, temp = x;
	while ( y )
	{
		if (y & 1)
			Res *= temp;
		temp *= temp;
		y /= 2;
	}
	return Res;
}
DB C[MAXN][MAXN];
void Prepare()
{
	C[0][0] = 1;
	for (Int i = 1; i < MAXN; ++ i)
	{
		C[i][0] = 1;
		for (Int j = 1; j <= i; ++ j)
			C[i][j] = C[i - 1][j - 1] + C[i - 1][j];
	}
	return ;
}
int main()
{
	int c, n, m;
	Prepare();
	while (~ scanf("%d %d %d", &c, &n, &m) && c)
	{
		if (m > c || m > n || (n - m) & 1)
		{
			printf("0.000\n");
			continue;
		}
		DB fuck = 0;
		for (Int i = 0; i <= m; ++ i)
			for (Int j = 0; j <= c - m; ++ j)
				fuck += ((i & 1) ? -1 : 1) * C[m][i] * C[c - m][j] * QuickPow((c - 2 * i - 2 * j) * 1.0 / c, n);
		fuck *= C[c][m];
		fuck /= QuickPow(2.0, c);
		printf("%.3lf\n", fuck);
	}
	return 0;
}