LOJ2047 气氛

气氛

有人向住在北大街的这 \(n\) 个人提了 \(n-1\) 个问题，根据每个人的回答，他会被分配一个 \(n-1\) 维的零一坐标，也就是一个点。这样 \(n\) 个点可以恰好构成一个 \(n-1\) 维空间中的凸包。

北大街的居民认为，在这个多面体内，便是华夏；多面体之外，便是蛮夷。我们可以很容易的计算出华夏部分的广义凸包体积。

有一天，清华路的 B 君来北大街玩，听说了这个故事觉得很有趣，于是也试着给出了这 \(n - 1\) 个问题的答案，

清华路的 B 君，当然认为自己属于华夏，但是北大街表示在 \(n-1\) 维空间中如果有 \(n+1\) 个点的话，华夏部分的体积难以计算。

这下子气氛突然江化。所以这个问题就留给你了，输入 \(n - 1\) 维度空间中的 \(n + 1\) 个点，求广义凸包的体积。

由于这个体积可能不是整数，你只需要输出体积乘以 \(n - 1\) 的阶乘，然后对 \(1000000007\) 取模的结果。

\(1 \leq t \leq 100\)，\(3 \leq n \leq 35\)，点的坐标一定是 \(0\) 或者 \(1\)。

题解

http://jklover.hs-blog.cf/2020/05/30/Loj-2407-气氛-Luft/#more

高维多面体的体积计算 + 高斯消元的一个 trick.

如果只有 \(n\) 个 \(n-1​\) 维的点,那么凸包其实就是它们围成的高维多面体.

它的体积就是选一个起点, 将 \(n-1\) 个 \(n-1\) 维向量形成的矩阵的行列式绝对值除掉 \((n-1)!​\) .

现在有 \(n+1\) 个点,求它们形成的凸包体积.

有一个通用的结论,这 \(n+1\) 个 \(n-1\) 维点的凸包体积等于任选出 \(n\) 个点的高维多面体体积之和 \(/2​\) .

枚举哪一个点没选,用高斯消元求出其余 \(n​\) 个点形成的高维多面体体积即可.

注意不能直接在模意义下计算,否则只能算出行列式在模意义下的值,但无法判定符号,得出绝对值.

所以必须要用真实值消元,为了避免 double 的精度误差,可以将真实值和模意义下的值分别维护.

最后根据真实值的符号决定是否对模意义下的值取相反数.

时间复杂度 \(O(t\cdot n^4)\) .

O int N=40;
int n,p[N][N],b[N][N];
float128 a[N][N];

int solve(){
	int ans=1;
	float128 res=1;
	for(int i=1;i<=n-1;++i){
		int p=i;
		for(int j=i;j<=n-1;++j)if(b[j][i]) {p=j; break;}
		if(!b[p][i]) return 0;
		if(p!=i){
			swap(b[p],b[i]),swap(a[p],a[i]);
			ans=mod-ans,res=-res;
		}
		ans=mul(ans,b[i][i]),res*=a[i][i];
		int inv=fpow(b[i][i],mod-2);
		for(int j=i+1;j<=n-1;++j)if(b[j][i]){
			int coef=mul(mod-b[j][i],inv);
			float128 tmp=-a[j][i]/a[i][i];
			for(int k=i;k<=n-1;++k)
				b[j][k]=add(b[j][k],mul(b[i][k],coef)),a[j][k]+=a[i][k]*tmp;
		}
	}
	if(res<0) ans=mod-ans;
	return ans;
}
void real_main(){
	read(n);
	for(int i=1;i<=n+1;++i)
		for(int j=1;j<=n-1;++j) read(p[i][j]);
	int ans=0;
	for(int ban=1;ban<=n+1;++ban){
		int st=ban==1?2:1;
		int tot=0;
		for(int i=1;i<=n+1;++i)if(i!=ban and i!=st){
			++tot;
			for(int j=1;j<=n-1;++j)
				a[tot][j]=p[i][j]-p[st][j],b[tot][j]=add(p[i][j],mod-p[st][j]);
		}
		ans=add(ans,solve());
	}
	ans=mul(ans,i2);
	printf("%d\n",ans);
}
int main(){
	for(int T=read<int>();T--;) real_main();
	return 0;
}