[Luogu P4035] [BZOJ 1013] [JSOI2008] 球形空间产生器

洛谷传送门

BZOJ传送门

题目描述

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

输入输出格式

输入格式：

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

输出格式：

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

输入输出样例

输入样例#1：

2
0.0 0.0
-1.0 1.0
1.0 0.0

输出样例#1：

0.500 1.500

说明

提示：给出两个定义：

1. 球心：到球面上任意一点距离都相等的点。
2. 距离：设两个$n$维空间上的点$A, B$的坐标为 $(a_1, a_2, \cdots , a_n), (b_1, b_2, \cdots , b_n)$，则$AB$的距离定义为： $dist = \sqrt{ (a_1-b_1)^2 + (a_2-b_2)^2 + \cdots + (a_n-b_n)^2 }$。

   解题分析

   高斯消元模板题…

   根据题意， 我们会有$N+1$个未知数， 即球心的坐标和半径。 我们可以列出形如以下的$N+1$个方程：
   

   $$R^2=(x_1-a_1)^2+(x_2-a_2)^2+...+(x_n-a_n)^2$$
   其中 $x_i$表示球心 $i$维的坐标。

   这玩意不好解， 我们拆开以一下：
   

   $$R^2-\sum_{i=1}^nx_i^2=\sum_{i=1}^{n}a_i^2-2\times \sum_{i=1}^{n}a_ix_i$$
   这下左边 $N+1$个方程作差可以得到0， 右边的 $\sum_{i=1}^na_i^2$为常数， 这下就顺利转变成了 N N <script type="math/tex" id="MathJax-Element-28">N</script>个一次方程组， 可以解出来了。

   代码如下：

   #include <cstdio>
   #include <cstring>
   #include <cctype>
   #include <cstdlib>
   #include <cmath>
   #include <algorithm>
   #define R register
   #define IN inline
   #define gc getchar()
   #define W while
   #define MX 15
   #define db double
   db dat[MX][MX], mat[MX][MX], sum[MX];
   int dim, bd;
   void Gauss()
   {
       R int i, j, k;
       for (i = 1; i <= dim; ++i)
       {
           for (j = i + 1; j <= dim; ++j)
           if(fabs(mat[j][i]) > fabs(mat[i][i])) std::swap(mat[i], mat[j]);
           for (j = i + 1; j <= bd; ++j) mat[i][j] /= mat[i][i];
           for (j = 1; j <= dim; ++j)
           {
               if(i ^ j)
               {
                   for (k = i + 1; k <= bd; ++k)
                   mat[j][k] -= mat[j][i] * mat[i][k];
               }
           }
       }
   }
   int main(void)
   {
       scanf("%d", &dim), bd = dim + 1;
       for (R int i = 0; i <= dim; ++i)
       for (R int j = 1; j <= dim; ++j)
       scanf("%lf", &dat[i][j]), sum[i] += dat[i][j] * dat[i][j];
       for (R int i = 1; i <= dim; ++i)
       {
           for (R int j = 1; j <= dim; ++j)
           mat[i][j] = 2 * (dat[i][j] - dat[i - 1][j]);
           mat[i][bd] = sum[i] - sum[i - 1];
       }
       Gauss();
       for (R int i = 1; i <= dim; ++i) printf("%.3lf ", mat[i][bd]);
   }