TopCoder 17403 See All Differences

这题和 https://atcoder.jp/contests/abc189/tasks/abc189_f 是相似的。

首先我们设/(f(mask,number)/)表示考虑若干个数，目前出现的差在/(mask/)中，最后一个数是/(number/)时，数列的期望长度。

但是我们发现这样比较难以转移，因为我们不知道转移过来的状态是否是合法的。（比如当你是rolled最后一个数的时候，mask只能包含rolled中出现的差，不能包含其它的数）

那么，我们不妨倒过来考虑——设/(f(mask,number)/)表示后面的若干个数中（因为数组是可以无限长的，故不存在什么合法不合法的问题），出现的差为mask，最前面的数为number时的答案。

则最终答案就是/(f(begin,rolled.back())/)，其中/(begin/)包含了所有不在/(rolled/)中的差。

于是我们就可以得到转移方程式：

/[f(mask,number)= /frac{1}{D} /sum_{number’=1}^{D} f(mask’,number’)+1
/]

若/(ch=|number-number’|/)在/(mask/)中，则/(mask’/)等于/(mask/)除去/(ch/)后的集合；否则/(mask’/)等于/(mask/)。（表示之前出现过了）

那么现在有一个问题——/(mask’/)是有可能等于/(mask/)的，那转移必将出现循环，该如何解决？

考虑对于每一个/(number/)，都会得到像上面转移那样的方程，故/(D/)个不同的/(number/)就构成了一个方程组。我们对其高斯消元，就可以解出每一个/(f(mask,number)/)了。

具体的，我们从小到大枚举/(mask/)（但/(f(0,number)=0/)是已知的，故/(mask/)要从/(1/)枚举），如果/(mask’=mask/)，那么将/(number’/)的系数减去/(/frac{1}{D}/)；否则将常数项加上/(/frac{1}{D} f(mask’,number’)/)，并且常数项一开始的值为/(1/)（因为转移后面有一个/(+1/)）。

时间复杂度为/(O(2^DD^3)/)。

#include<bits/stdc++.h>
#define debug(...) std::cerr<<#__VA_ARGS__<<" : "<<__VA_ARGS__<<std::endl

int pos[18],vis[18];
double a[18][18],f[1<<18][18];

void gauss(int n) {
	memset(pos,0,sizeof pos);
	memset(vis,0,sizeof vis);
	for(int i=1;i<=n;i++) {
		int maxp=0;
		for(int j=1;j<=n;j++)
			if((!maxp||abs(a[maxp][i])<abs(a[j][i]))&&!vis[j])
				maxp=j;
		if(a[maxp][i]==0) {
			puts("No Solution");
			exit(0);
		}
		pos[i]=maxp,vis[maxp]=1;
		double tim=a[maxp][i];
		for(int j=i;j<=n+1;j++)
			a[maxp][j]/=tim;
		for(int j=1;j<=n;j++) {
			if(j==maxp) continue;
			double tim=a[j][i];
			for(int k=i;k<=n+1;k++)
				a[j][k]-=a[maxp][k]*tim;
		}
	}
}

struct SeeAllDifferences {
	double solve(int D,std::vector<int> rolled) {
		for(int mask=1;mask<(1<<D);mask++) {
			memset(a,0,sizeof a);
			for(int x=0;x<D;x++) a[x+1][x+1]=a[x+1][D+1]=1;
			for(int x=0;x<D;x++) for(int y=0;y<D;y++) {
				int diff=abs(x-y);
				if(mask&(1<<diff)) {
					a[x+1][D+1]+=f[mask^(1<<diff)][y]/(double)D;
				} else {
					a[x+1][y+1]-=1.0/(double)D;
				}
			}
			gauss(D);
			for(int x=0;x<D;x++) f[mask][x]=a[pos[x+1]][D+1];
		}
		int mask=(1<<D)-1;
		for(int i=1;i<(int)rolled.size();i++) {
			int diff=abs(rolled[i]-rolled[i-1]);
			if(mask&(1<<diff)) mask^=1<<diff;
		}
		return f[mask][rolled.back()-1];
	}
};