BZOJ 1013: [JSOI2008]球形空间产生器sphere 高斯消元

最新推荐文章于 2020-02-12 16:41:03 发布

原创最新推荐文章于 2020-02-12 16:41:03 发布 · 558 阅读

0 ·

CC 4.0 BY-SA版权

bzoj 同时被 2 个专栏收录

88 篇文章

订阅专栏

数学

23 篇文章

订阅专栏

本文介绍了一种使用高斯消元法求解n维空间中通过n+1个点确定的球心坐标的方法。通过将n+1个点间距离方程转化为n个线性方程组，进而利用高斯消元法求解未知的球心坐标。

不停地把第i列上的未知数消去。

* 裸的高斯消元题
* 二维平面上的圆上的点与圆心的距离有(x-a)^2+(y-b)^2 = r^2
* 三维空间上的球上的点与球心的距离有(x-a)^2+(y-b)^2+(z-c)^2 = r^2
* 同理：在n维空间上的球上的点与球心的距离有sigma((xi-ai)^2) = r^2，圆心为(a1,a2,...,an)
*
* 另外，在二维平面上，可由三点（不共线）确定一个园，在三维上四点（不共线）确定一个球，
* 同理，在n维平面上，可由n+1个点（不共线）确定一个n维的球。
*
* 这样，题目就可以转化为n+1个方程组，但是是平方级别的，如何转化为一维的？
* 我们不妨对于相邻的两个方程组左右分别相减，可以发现：
* 2*(x21 - x11)*x1 + 2*(x22 - x12)*x2 +...+2*(x2n - x1n) = (x21^2 - x11^1)+...+(x2n^2 - x1n^2)
* 这样，由n+1个方程组就可以转化为了n个一维的方程组了。下面，直接用高斯消元法即可解决该问题

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
#define eps 1e-8
int n,m;
double save[100][100],a[100][100],x[100];
void guass()
{
    m=n+1;
    int row,col,i,j,k;
    for(row=1,col=1;row<n,col<m;row++,col++)
    {
        k=row;
        for(int i=row+1;i<=n;i++) if(fabs(a[i][col])>fabs(a[k][col])) k=i;

        if(k!=row)
        {
            for(i=col; i<=m; i++)  swap(a[k][i],a[row][i]);
        }
        for(i=row+1; i<=n; i++)
        {
            if(fabs(a[i][col])<eps)
                continue;
            double t=a[i][col]/a[row][col];
            a[i][col]=0.0;
            for(int j=col+1;j<=m;j++)
                a[i][j]-=a[row][j]*t;
        }
    }
    for(i=n;i>=1;i--)
    {
        x[i]=a[i][m];
        for(j=i+1;j<=n;j++)
            x[i]-=x[j]*a[i][j];
        x[i]/=a[i][i];
    }
}
int main()
{
    scanf("%d",&n);
    for(int i=1;i<=n+1;i++)
    {
        for(int j=1;j<=n;j++)
        {
            scanf("%lf",&save[i][j]);
        }
    }
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<=n;j++)
        {
            a[i][j]=2.0*(save[i][j]-save[i+1][j]);
            a[i][n+1]+=save[i][j]*save[i][j]-save[i+1][j]*save[i+1][j];
        }
    }
    guass();
    for(int i=1;i<n;i++)
    {
        printf("%.3f ",x[i]);
    }
    printf("%.3f\n",x[n]);
    return 0;
}