[bzoj 1013] [JSOI2008]球形空间产生器sphere：高斯消元

题意：n维空间中有一个球面，给(n+1)个球面上的点的坐标，求球心坐标。n维空间中的球面定义为到定点（球心）的欧几里德距离等于定值（半径）的点的集合。数据保证有解。

看起来可以待定系数法，但是会不会产生二次项呢？

对于每个点(a1,a2,…,an)，可以列方程：

(x1−a1)2+(x2−a2)2+⋯+(xn−an)2=r2

把含未知量的项放到一边，常数项放到另一边，得：

2a1x1+2a2x2+⋯2anxn−(x21+x22+⋯+x2n)+r2=a21+a2+⋯+a2n

如果把括号内的一串平方和看作整体，也还是有(n+2)个变量，(n+1)个方程搞不定诶……

但是每个方程都含有−(x21+x22+⋯+x2n)+r2，一作差不就消掉了吗？这样有n个方程，n个未知数，不错。

等等……把−(x21+x22+⋯+x2n)+r2看作一个整体不就好了吗？

整体思想掌握地比我好的读者请略过上述三段。

其实都可以啦……黄学长是作的差。这样可以减少一个未知数。

为了方便行的交换，我把增广矩阵的每一行申请为动态内存。

#include <cstdio>
#include <algorithm>
#include <cmath>
using namespace std;
const double eps = 1e-9;
const int MAX_N = 10;
double* M[MAX_N+1];
int n;

void gauss()
{
for (int i = 0, j; i <= n; ++i) {
for (j = i; j <= n && fabs(M[j][i]) < eps; ++j)
;
swap(M[i], M[j]);
double t = M[i][i];
M[i][i] = 1;
for (j = i+1; j <= n+1; ++j)
M[i][j] /= t;
for (j = 0; j <= n; ++j) {
if (j == i)
continue;
double r = M[j][i];
M[j][i] = 0;
for (int k = i+1; k <= n+1; ++k)
M[j][k] -= r*M[i][k];
}
}
}

int main()
{
scanf("%d", &n);
for (int i = 0; i <= n; ++i) {
M[i] = new double[n+2]();
M[i][0] = 1;
for (int j = 1; j <= n; ++j) {
scanf("%lf", &M[i][j]);
M[i][n+1] += M[i][j] * M[i][j];
M[i][j] *= 2;
}
}
gauss();
for (int i = 1; i <= n; ++i)
printf("%.3f%c", M[i][n+1], " \n"[i == n]);
return 0;
}