Skip to content

BZOJ - 1013. 球形空间产生器

检测到 KaTeX 加载失败,可能会导致文中的数学公式无法正常渲染。

#题面

#题目描述

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

#输入格式

第一行是一个整数 nn。接下来的 n+1n + 1 行,每行有 nn 个实数,表示球面上一点的 nn 维坐标。每一个实数精确到小数点后 66 位,且其绝对值都不超过 2000020000

#输出格式

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

#输入输出样例

样例输入 #1

2
0.0 0.0
-1.0 1.0
1.0 0.0

样例输出 #1

0.500 1.500

#数据范围与约定

对于 100%100\% 的数据,1n101 \leq n \leq 10

#提示

给出两个定义:

  1. 球心:到球面上任意一点距离都相等的点。
  2. 距离:设两个 n 维空间上的点 A,BA, B 的坐标为 (a1,a2, ,an),(b1,b2, ,bn)(a_1, a_2, \cdots , a_n), (b_1, b_2, \cdots , b_n),则 ABAB 之间的距离定义为:dist=i=1n(aibi)2\mathit{dist} = \sqrt{ \sum_{i = 1}^{n} (a_i - b_i)^2 }

#思路

设球的半径为 RR,球心坐标为 (x1,x2,,xn)(x_1, x2, \ldots, x_n)。那么根据球体的性质,显然有:

j=0n(ai,jxj)2=R2\sum_{j = 0}^{n} (a_{i, j} - x_j)^2 = R^2

那么展开来写可以得到一个方程组的形式:

{(a1,1x1)2+(a1,2x2)2++(a1,nxn)2=R2(a2,1x1)2+(a2,2x2)2++(a2,nxn)2=R2(an+1,1x1)2+(an+1,2x2)2++(an+1,nxn)2=R2\begin{cases} \begin{aligned} (a_{1, 1} - x_1)^2 + (a_{1, 2} - x_2)^2 + \cdots + (a_{1, n} - x_n)^2 &= R^2 \\ (a_{2, 1} - x_1)^2 + (a_{2, 2} - x_2)^2 + \cdots + (a_{2, n} - x_n)^2 &= R^2 \\ \vdots \\ (a_{n + 1, 1} - x_1)^2 + (a_{n + 1, 2} - x_2)^2 + \cdots + (a_{n + 1, n} - x_n)^2 &= R^2 \\ \end{aligned} \end{cases}

但这个方程组是由 n+1n + 1nn 元二次方程构成的,不是线性方程组。本着「出题人绝对不会多给 1 Byte 有用数据」的原则,我们可以察觉到 n+1n + 1 中的异样之处。通过将相邻两个方程作差, 得到的结果是:

{2(a2,1a1,1)x1+2(a2,2a1,2)x2++2(a2,na1,n)xn=a2,12+a2,22++a2,n2a1,12a1,22a1,n22(a3,1a2,1)x1+2(a3,2a2,2)x2++2(a3,na2,n)xn=a3,12+a3,22++a3,n2a2,12a2,22a2,n22(an+1,1an,1)x1+2(an+1,2an,2)x2++2(an+1,nan,n)xn=an+1,12+an+1,22++an+1,n2an,12an,22an,n2\begin{cases} \begin{aligned} 2(a_{2, 1} - a_{1, 1}) x_1 + 2(a_{2, 2} - a_{1, 2}) x_2 + \cdots + 2(a_{2, n} - a_{1, n}) x_n &= a_{2, 1}^2 + a_{2, 2}^2 + \cdots + a_{2, n}^2 - a_{1, 1}^2 - a_{1, 2}^2 - \cdots - a_{1, n}^2 \\ 2(a_{3, 1} - a_{2, 1}) x_1 + 2(a_{3, 2} - a_{2, 2}) x_2 + \cdots + 2(a_{3, n} - a_{2, n}) x_n &= a_{3, 1}^2 + a_{3, 2}^2 + \cdots + a_{3, n}^2 - a_{2, 1}^2 - a_{2, 2}^2 - \cdots - a_{2, n}^2 \\ \vdots \\ 2(a_{n + 1, 1} - a_{n, 1}) x_1 + 2(a_{n + 1, 2} - a_{n, 2}) x_2 + \cdots + 2(a_{n + 1, n} - a_{n, n}) x_n &= a_{n + 1, 1}^2 + a_{n + 1, 2}^2 + \cdots + a_{n + 1, n}^2 - a_{n, 1}^2 - a_{n, 2}^2 - \cdots - a_{n, n}^2 \\ \end{aligned} \end{cases}

可以转化成增广矩阵的形式:

[2(a2,1a1,1)2(a2,2a1,2)2(a2,na1,n)j=1n(a2,j2a1,j2)2(a3,1a2,1)2(a3,2a2,2)2(a3,na2,n)j=1n(a3,j2a2,j2)2(an+1,1an,1)2(an+1,2an,2)2(an+1,nan,n)j=1n(an+1,j2an,j2)]\left[ \begin{array}{cccc|c} 2(a_{2, 1} - a_{1, 1}) & 2(a_{2, 2} - a_{1, 2}) & \cdots & 2(a_{2, n} - a_{1, n}) & \sum_{j = 1}^{n} (a_{2, j}^2 - a_{1, j}^2) \\ 2(a_{3, 1} - a_{2, 1}) & 2(a_{3, 2} - a_{2, 2}) & \cdots & 2(a_{3, n} - a_{2, n}) & \sum_{j = 1}^{n} (a_{3, j}^2 - a_{2, j}^2) \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 2(a_{n + 1, 1} - a_{n, 1}) & 2(a_{n + 1, 2} - a_{n, 2}) & \cdots & 2(a_{n + 1, n} - a_{n, n}) & \sum_{j = 1}^{n} (a_{n + 1, j}^2 - a_{n, j}^2) \\ \end{array} \right]

对该矩阵进行高斯消元即可。

#代码

C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
#include <algorithm>
#include <cmath>
#include <iomanip>
#include <iostream>

using std::cin;
using std::cout;
const char endl = '\n';

const int N = 15;
const double eps = 1e-6;

int n;
double a[N][N], b[N][N];

int main() {
std::ios::sync_with_stdio(false);

cin >> n;

for (int i = 1; i <= n + 1; i++) {
for (int j = 1; j <= n; j++) {
cin >> a[i][j];
}
}

for (int i = 1; i <= n; i++) {
for (int j = 1; j <= n; j++) {
b[i][j] = 2.0 * (a[i][j] - a[i + 1][j]);
b[i][n + 1] += std::pow(a[i][j], 2) - std::pow(a[i + 1][j], 2);
}
}

for (int i = 1; i <= n; i++) {
for (int j = i; j <= n; j++) {
if (std::abs(b[j][i]) > eps) {
std::swap(b[i], b[j]);
}
}

for (int j = 1; j <= n; j++) {
if (i == j) continue;
double x = b[j][i] / b[i][i];
for (int k = i; k <= n + 1; k++) {
b[j][k] -= b[i][k] * x;
}
}
}

for (int i = 1; i <= n; i++) {
cout << std::fixed << std::setprecision(3) << b[i][n + 1] / b[i][i] << ' ';
}
cout << endl;

return 0;
}