$$\Large\texttt{My Blog}$$


Description

题目链接:Luogu 4035

有一个球形空间产生器能够在 $n$ 维空间中产生一个坚硬的球体。现在,你被困在了这个 $n$ 维球体中,你只知道球面上 $n+1$ 个点的坐标 $(x_{i,1},x_{i,2},\dots,x_{i,n})$ ,你需要以最快的速度确定这个 $n$ 维球体的球心坐标,以便于摧毁这个球形空间产生器。答案保留 $3$ 位小数。

数据范围: $1\le n\le 10$ , $\vert x_{i,j}\vert\le 20000$


Solution

得到 $n+1$ 个等式: $$\begin{cases}(x_1-p_{1,1})^2+(x_2-p_{1,2})^2+(x_3-p_{1,3})^2+\cdots+(x_n-p_{1,n})^2=R^2 \\(x_1-p_{2,1})^2+(x_2-p_{2,2})^2+(x_3-p_{2,3})^2+\cdots+(x_n-p_{2,n})^2=R^2 \\\vdots \\(x_1-p_{n,1})^2+(x_2-p_{n,2})^2+(x_3-p_{n,3})^2+\cdots+(x_n-p_{n,n})^2=R^2 \\(x_1-p_{n+1,1})^2+(x_2-p_{n+1,2})^2+(x_3-p_{n+1,3})^2+\cdots+(x_n-p_{n+1,n})^2=R^2\end{cases}$$ 然后用第 $n+1$ 个等式去消去前 $n$ 个等式的 ${x_i}^2$ 项,这样就得到了一个标准的 $n$ 元 $1$ 次方程组,直接高斯消元求出 $n$ 个未知数。

时间复杂度: $O(n^3)$


Code

#include <cstdio>
#include <cmath>
#include <algorithm>

const int N=15;
int n;
double p[N][N],a[N][N],b[N],x[N];

void Gauss(int n) {
    for(int i=1;i<=n;++i) {
        int p=i;
        for(int k=i+1;k<=n;++k) if(fabs(a[k][i])>fabs(a[p][i])) p=k;
        if(i!=p) std::swap(a[i],a[p]),std::swap(b[i],b[p]);
        for(int k=i+1;k<=n;++k) {
            double d=a[k][i]/a[i][i];
            b[k]-=d*b[i];
            for(int j=i;j<=n;++j) a[k][j]-=d*a[i][j];
        }
    }
    for(int i=n;i>=1;--i) {
        for(int j=i+1;j<=n;++j) b[i]-=x[j]*a[i][j];
        x[i]=b[i]/a[i][i];
    }
}
void init() {
    for(int i=1;i<=n;++i) for(int j=1;j<=n;++j) {
        a[i][j]=2.0*(p[n+1][j]-p[i][j]),b[i]+=p[n+1][j]*p[n+1][j]-p[i][j]*p[i][j];
    }
}
int main() {
    scanf("%d",&n);
    for(int i=1;i<=n+1;++i) for(int j=1;j<=n;++j) scanf("%lf",&p[i][j]);
    init();
    Gauss(n);
    for(int i=1;i<=n;++i) printf("%.3lf%c",x[i]," \n"[i==n]);
    return 0;
}