[bzoj 1013][jsoi 2008] 球形空间产生器

题目大意

给定 (n) 维空间上一个球面上的 (n + 1) 个点,求其圆心。数据保证有解。

(1 leqslant n leqslant 10)

(|坐标| leqslant 20,000)

题目链接

BZOJ 1013

题解

考虑两个点 (P)(Q) ,记圆心为 (O) ,有 [
begin{align}
sum_{i = 1}^{n} (P_i - Oi)^2 &= sum{i = 1}^{n} (Q_i - Oi)^2 \
sum
{i = 1}^{n} P_i^2 - 2 P_i Oi &= sum{i = 1}^{n} Q_i^2 - 2 Q_i Oi \
sum
{i = 1}^{n} 2(P_i - Q_i)Oi &= sum{i = 1}^{n} P_i^2 - Q_i^2
end{align}
]
发现二次项被约掉,于是可以高斯消元解方程。

代码

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

#include <cmath>
#include <algorithm>
const int MAXN = 10;
const double EPS = 1e-7;
double a[MAXN][MAXN + 1];
bool (int n) {
for (int i = 0; i < n; i++) {
int max = i;
for (int j = i + 1; j < n; j++) if (fabs(a[j][i]) > fabs(a[max][i])) max = j;
if (fabs(a[max][i]) < EPS) return false;
if (max != i) for (int j = i; j <= n; j++) std::swap(a[i][j], a[max][j]);
for (int j = 0; j < n; j++) if (i != j) {
for (int k = n; k >= i; k--) a[j][k] -= a[i][k] / a[i][i] * a[j][i];
}
}
return true;
}
int main() {
int n;
scanf("%d", &n);
for (int i = 0; i <= n; i++) for (int j = 0; j < n; j++) {
double x;
scanf("%lf", &x);
if (i != n) a[i][j] += 2 * x, a[i][n] += x * x;
if (i != 0) a[i - 1][j] -= 2 * x, a[i - 1][n] -= x * x;
}
if (!GaussJordan(n)) return puts("-1"), 0;
for (int i = 0; i < n; i++) printf("%.3lf%c", a[i][n] / a[i][i], i == n - 1 ? 'n' : ' ');
return 0;
}