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

Description

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

Input & Output

Input

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

Output

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

Sample

懒得粘了

Solution

空间距离地计算其实给了我们一些灵感，不妨设一个球心，让所有店到球心距离相等，但是这样构造不出矩阵的最右常数列，而给出的第n+1个坐标提供了帮助。把根号去掉取等，然后胡乱展开再移项，高斯消元即可。具体等式看代码就理解了……蒟蒻不会码公式……
Code：

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>

using std::min;
using std::max;
using std::swap;
using std::fabs;

const int maxn = 15;

double data[maxn][maxn],matrix[maxn][maxn];
int n;

double pow(double x)
{
    return x * x;
}
void Gauss()
{
    int cur;
    for(int i = 1; i <= n; ++i)
    {
        cur = i;
        for(int j = i + 1; j <= n; ++j)
            fabs(matrix[cur][i]) < fabs(matrix[j][i]) ? cur = j : 1;
        for(int j = i; j <= n + 1; ++j) swap(matrix[cur][j],matrix[i][j]);
        for(int j = i + 1; j <= n + 1; ++j)
            matrix[i][j] /= matrix[i][i];
        matrix[i][i] = 1;
        for(int j = i + 1; j <= n; ++j)
        {
            for(int k = i + 1; k <= n + 1; ++k)
                matrix[j][k] -= matrix[j][i] * matrix[i][k];
            matrix[j][i] = 0;
        }
    }
    for(int i = n; i >= 1; --i)
    {
        for(int j = i + 1; j <= n; ++j)
            matrix[i][n+1] -= matrix[i][j] * matrix[j][n+1], matrix[i][j] = 0;
        matrix[i][n+1] /= 1;
        matrix[i][i] = 1;
    }
    for(int i = 1; i <= n; ++i)
        printf("%.3lf ",matrix[i][n+1]);
    putchar('\n');
}
void init()
{
    for(int i = 1; i <= n; ++i)
    {
        for(int j = 1; j <= n; ++j)
            matrix[i][j] = -2 * data[i][j] + 2 * data[n+1][j];
        for(int j = 1; j <= n; ++j)
            matrix[i][n+1] += pow(data[n+1][j]) - pow(data[i][j]);
    }
}

int main()
{
    scanf("%d",&n);
    for(int i = 1; i <= n + 1; ++i)
        for(int j = 1; j <= n; ++j)
            scanf("%lf",&data[i][j]);
    init();
    Gauss();
    return 0;
}