##[$AcWing$ $207$. 球形空间产生器 ](https://www.acwing.com/problem/content/description/209/) ### 一、题目描述 有一个球形空间产生器能够在 $n$ 维空间中产生一个坚硬的球体。 现在,你被困在了这个 $n$ 维球体中,你只知道球面上 $n+1$ 个点的坐标,你需要以最快的速度确定这个 $n$ 维球体的球心坐标,以便于摧毁这个球形空间产生器。 **注意**: 数据保证有唯一解。 **输入格式** 第一行是一个整数 $n$。 接下来的 $n+1$ 行,每行有 $n$ 个实数,表示球面上一点的 $n$ 维坐标。 每一个实数精确到小数点后 $6$ 位,且其绝对值都不超过 $20000$。 **输出格式** 有且只有一行,依次给出球心的 $n$ 维坐标($n$ 个实数),两个实数之间用一个空格隔开。 每个实数精确到小数点后 $3$ 位。 **数据范围** $1≤n≤10$ **输入样例:** ```cpp {.line-numbers} 2 0.0 0.0 -1.0 1.0 1.0 0.0 ``` **输出样例:** ```cpp {.line-numbers} 0.500 1.500 ``` ### 二、解题思路 **[前置知识:高斯消元法求解线性方程组](https://www.cnblogs.com/littlehb/p/15385357.html)** 设球心坐标为$(x_1,x_2,...,x_n)$ 对于给定的$n+1$个点坐标为$(b_{i1},b_{i2},...,b_{in})$ 根据题意有如下$n+1$个方程构成的方程组 $$ \large \left\{\begin{matrix} (b_{01}-x_1)^2 +(b_{02}-x_2)^2 +...+(b_{0n}-x_n)^2=R^2& \ ① \\ (b_{11}-x_1)^2 +(b_{12}-x_2)^2 +...+(b_{1n}-x_n)^2=R^2& \ ② \\ (b_{21}-x_1)^2 +(b_{22}-x_2)^2 +...+(b_{2n}-x_n)^2=R^2& \ ③ \\ ... \\ (b_{n1}-x_1)^2 +(b_{n2}-x_2)^2 +...+(b_{nn}-x_n)^2=R^2& \ \end{matrix}\right. $$ 因为我们学过的知识,只有 **高斯消元可以解这样的方程组**,但本题的方程是二次方的,肯定是不行,**需要想办法把二次方消掉**,化简为一次的多元线性方程组,才能使用高斯消元。 进行公式变换,把二项全部消掉: $\large \displaystyle \left\{\begin{matrix} ②-① \\ ③-① \\ ... \\ (n+1) -① & \end{matrix}\right.$ 有: $$ \large \displaystyle \left\{\begin{matrix} 2(b_{11}-b_{01})x_1+2(b_{12}-b_{02})x_2+...+2(b_{1n}-b_{0n})x_n=b_{11}^2+b_{12}^2+...+b_{1n}^2-b_{01}^2-b_{02}^2-...b_{0n}^2 \\ 2(b_{21}-b_{01})x_1+2(b_{22}-b_{02})x_2+...+2(b_{2n}-b_{0n})x_n=b_{11}^2+b_{12}^2+...+b_{1n}^2-b_{01}^2-b_{02}^2-...b_{0n}^2 \\ ... \\ 2(b_{n1}-b_{01})x_1+2(b_{n2}-b_{02})x_2+...+2(b_{n1}-b_{0n})x_n=b_{n1}^2+b_{n2}^2+...+b_{nn}^2-b_{01}^2-b_{02}^2-...b_{0n}^2 \\ \end{matrix}\right. $$ 令: $$\large \displaystyle \left\{\begin{matrix} 2(b_{ij}-b_{0j})=a_{ij} \\ b_{i1}^2+b_{i2}^2+...+b_{in}^2-b_{01}^2-b_{02}^2-...b_{0n}^2=a_{i,{n+1}} \end{matrix}\right. \ i,j \in [1,n] $$ 则有: $$ \large \left\{\begin{matrix} a_{11}x_1 +a_{12}x_2+...+a_{1n}x_n=a_{1,{n+1}} \\ a_{21}x_1 +a_{22}x_2+...+a_{2n}x_n=a_{2,{n+1}} \\ ... \\ a_{n1}x_1 +a_{n2}x_2+...+a_{nn}x_n=a_{n,{n+1}} \end{matrix}\right. $$ 如此变化后,就得到一个$n$个变量,$n$个等式的 **线程方程组**,可以使用 **高斯消元** 来处理。 ### 三、实现代码 ```cpp {.line-numbers} #include using namespace std; const int N = 15; const double eps = 1e-8; // 小数的精度值 int n; double a[N][N], b[N][N]; int gauss() { // 高斯消元,答案存于a[i][n]中,0 <= i < n int r = 1; // 先按行后按列进行计算,当前行是第1行 for (int c = 1; c <= n; c++) { // 枚举每一列 int t = r; // 防破坏r,复制出t for (int i = r; i <= n; i++) // 当前行需要找它的后续行 if (abs(a[i][c]) > abs(a[t][c])) t = i; // t的任务是找出c列中系数最大值是哪一行 if (abs(a[t][c]) < eps) continue; // 如果c列绝对值最大的系数是0, 那么处理下一列 for (int i = c; i <= n + 1; i++) swap(a[t][i], a[r][i]); // 将绝对值最大的行与当前行交换 for (int i = n + 1; i >= c; i--) a[r][i] /= a[r][c]; // a[r][c]:行首系数,将当前行的行首通过除法变为1,倒序 for (int i = r + 1; i <= n; i++) // 用当前行r的c列,通过减法将后续行c列消成0 for (int j = n + 1; j >= c; j--) // 倒序,需要保留行首,逻辑和上面是一样的,行首值是变更系数,如果正序就把系数变成1了,后面就不对了 a[i][j] -= a[r][j] * a[i][c]; // a[i][c]:需要变化的乘法系数,减法:对位相消 r++; // 下一行 } if (r <= n) { // 如果没有成功执行完所有行,意味着中间存在continue,也就是某一列的系数都是0 for (int i = r; i <= n; i++) if (abs(a[i][n + 1]) > eps) return 0; // 系数是0,但结果不是0,无解 return 2; // 系数是0,结果也是0,x取啥都对,有无穷多组解 } // 代回求每个变量值 for (int i = n - 1; i; i--) // 行,倒序 for (int j = i + 1; j <= n; j++) // 列,倒三角,右上角应该都是0,对角线全是1 a[i][n + 1] -= a[i][j] * a[j][n + 1]; // 系数消为0 return 1; // 有唯一解 } int main() { cin >> n; for (int i = 0; i <= n; i++) // n+1行,下标从0开始 for (int j = 1; j <= n; j++) // n列 cin >> b[i][j]; // 转换为线性方程组 for (int i = 1; i <= n; i++) for (int j = 1; j <= n; j++) { a[i][j] += 2 * (b[i][j] - b[0][j]); a[i][n + 1] += b[i][j] * b[i][j] - b[0][j] * b[0][j]; } // Gauss消元 gauss(); // 输出答案 for (int i = 1; i <= n; i++) printf("%.3lf ", a[i][n + 1]); return 0; } ```