[JSOI2008]球形空间产生器sphere

Description

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

Input

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

Output

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

Sample Input

2
0.0 0.0
-1.0 1.0
1.0 0.0

Sample Output

0.500 1.500

HINT

提示：
给出两个定义：
1、 球心：到球面上任意一点距离都相等的点。
2、 距离：设两个n为空间上的点A,B的坐标为(a1,a2,…,an),(b1,b2,…,bn)，则AB的距离定义为：dist=(a1−b1)2+(a2−b2)2+…+(an−bn)2−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−√

Source

思路

高斯消元模板。
设球心O的坐标为(x1,x2,⋯,xn)，点P在j方向上的坐标为aP,j，则有：

|OA|=(aA,1−x1)2+(aA,2−x2)2+(aA,3−x3)2+⋯+(aA,n−xn)2=a2A,1−2aA,1x1+x21+⋯+a2A,n+2aA,nx2n+x2n

同理，

|OB|=a2B,1−2aB,1x1+x21+⋯+a2B,n+2aB,nx2n+x2n

又因为

|OA|=|OB|

所以

a2A,1−2aA,1x1+x21+⋯+a2A,n+2aA,nx2n+x2n=a2B,1−2aB,1x1+x21+⋯+a2B,n+2aB,nx2n+x2n

整理得

∑i=1n(2aA,i−2aB,i)xi=∑i=1n(a2A,i−a2B,i)

得出 n个方程，高斯消元解方程即可得出答案。

代码

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

const int maxn=10;
const double eps=0.0001;

int n;

struct matrix
{
  double a[maxn+1][maxn+2];

  inline int gauss()
  {
    for(register int  i=1; i<n; ++i)
      {
        int t=i;
        for(register int j=i+1; j<=n; ++j)
          {
            if(fabs(a[t][i])<fabs(a[j][i]))
              {
                t=j;
              }
          }
        if(fabs(a[t][i])<eps)
          {
            return 1;
          }
        for(register int j=i; j<=n+1; ++j)
          {
            std::swap(a[t][j],a[i][j]);
          }
        for(register int j=i+1; j<=n; ++j)
          {
            if(!a[j][i])
              {
                continue;
              }
            double f=a[j][i]/a[i][i];
            for(register int k=i; k<=n+1; ++k)
              {
                a[j][k]-=f*a[i][k];
              }
          }
      }
    for(register int i=n; i; --i)
      {
        a[i][n+1]=a[i][n+1]/a[i][i];
        a[i][i]=1;
        for(register int j=1; j<i; ++j)
          {
            a[j][n+1]-=a[j][i]*a[i][n+1];
          }
      }
    return 0;
  }

  inline int print_ans()
  {
    for(register int i=1; i<n; ++i)
      {
        printf("%.3lf ",a[i][n+1]);
      }
    printf("%.3lf",a[n][n+1]);
    return 0;
  }
};

matrix t;
double point[maxn+2][maxn+1];

int main()
{
  scanf("%d",&n);
  for(register int i=1; i<=n+1; ++i)
    {
      for(register int j=1; j<=n; ++j)
        {
          scanf("%lf",&point[i][j]);
        }
    }
  for(register int i=1; i<=n; ++i)
    {
      double sum=0;
      for(register int j=1; j<=n; ++j)
        {
          t.a[i][j]=2*(point[i][j]-point[i+1][j]);
          sum+=point[i][j]*point[i][j]-point[i+1][j]*point[i+1][j];
        }
      t.a[i][n+1]=sum;
    }
  t.gauss();
  t.print_ans();
  return 0;
}

转载于:https://www.cnblogs.com/Canopus-wym/p/10376247.html