[JSOI2008]球形空间产生器——[高斯消元]

【题目描述】

有一个球形空间产生器能够在 n 维空间中产生一个坚硬的球体。现在，你被困在了这个 n 维球体中，你只知道球面上 n+1 个点的坐标，你需要以最快的速度确定这个 n 维球体的球心坐标，以便于摧毁这个球形空间产生器。

【输入格式】

第一行是一个整数 n (1<=N=10)。接下来的 n+1 行，每行有 n 个实数，表示球面上一点的 n 维坐标。每一个实数精确到小数点后 6 位，且其绝对值都不超过 20000。

【输出格式】

有且只有一行，依次给出球心的 n 维坐标（ n 个实数），两个实数之间用一个空格隔开。每个实数精确到小数点后 3 位。数据保证有解。你的答案必须和标准输出一模一样才能够得分。

$SampleInput$

2
0.0 0.0
-1.0 1.0
1.0 0.0

$SampleOutput$

0.500 1.500

【题意分析】

给出$n+1$个有$n$维坐标的点，这些点到某个点的$n$维空间距离都相等，求这个点坐标。

那么设这个点坐标是$(x_1,x_2,...,x_n)$
第$1-n+1$个点坐标是$(a_{1,1},a_{1,2},...,a_{1,n}),(a_{2,1},a_{2,2},...,a_{2,n}),...,(a_{n+1,1},a_{n+1,2},...,a_{n+1,n}),$

由于这$n+1$个点到未知点的距离都相等

先设$n=2$，距离为$dist$

$\{\begin{array}{l}(a_{1,1}-x_1{)}^2+(a_{1,2}-x_2{)}^2=dis{t}^2\\ (a_{2,1}-x_1{)}^2+(a_{2,2}-x_2{)}^2=dis{t}^2\\ (a_{3,1}-x_1{)}^2+(a_{3,2}-x_2{)}^2=dis{t}^2\end{array}$

$\Rightarrow$ $\{\begin{array}{l}(a_{1,1}-x_1{)}^2+(a_{1,2}-x_2{)}^2=(a_{2,1}-x_1{)}^2+(a_{2,2}-x_2{)}^2\\ (a_{2,1}-x_1{)}^2+(a_{2,2}-x_2{)}^2=(a_{3,1}-x_1{)}^2+(a_{3,2}-x_2{)}^2\end{array}$

展开一下

$\Rightarrow$ $\{\begin{array}{l}2(a_{2,1}-a_{1,1})x1+2(a_{2,2}-a_{1,2})x2=a_{2,1}^2-a_{1,1}^2+a_{2,2}^2-a_{1,2}^2\\ 2(a_{3,1}-a_{2,1})x1+2(a_{3,2}-a_{2,2})x2=a_{3,1}^2-a_{2,1}^2+a_{3,2}^2-a_{2,2}^2\end{array}$

然后类比$n$的情况，高斯消元即可

Code:

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <cmath>
#include <algorithm>
#define MAXN 2000
using namespace std;

double data[MAXN][MAXN], mat[MAXN][MAXN], ans[MAXN];

int main () {
	int n; scanf ("%d", &n);
	for (register int i = 1; i <= n + 1; i++)
		for (register int j = 1; j <= n; j++)
			scanf ("%lf", &data[i][j]);
	for (register int i = 1; i <= n; i++)
		for (register int j = 1; j <= n; j++)
			mat[i][j] = 2 * (data[i + 1][j] - data[i][j]),
			mat[i][n + 1] += data[i + 1][j] * data[i + 1][j] - data[i][j] * data[i][j];
	for (register int i = 1; i <= n; i++) {
		int pos = i;
		for (register int j = i + 1; j <= n; j++)
			if (fabs (mat[j][i]) > fabs (mat[pos][i])) pos = j;
		swap (mat[i], mat[pos]);
		for (register int j = i + 1; j <= n; j++) {
			double alpha = mat[j][i] / mat[i][i];
			for (register int k = i; k <= n + 1; k++)
				mat[j][k] -= mat[i][k] * alpha;
		}
	}
	for (register int i = n; i >= 1; i--) {
		ans[i] = mat[i][n + 1];
		for (register int j = n; j >= 1; j--)
			if (j != i) ans[i] -= mat[i][j] * ans[j];
		ans[i] /= mat[i][i];
	}
	for (register int i = 1; i < n; i++) printf ("%.3lf ", ans[i]);
	printf ("%.3lf", ans[n]);
	return 0;
}