[bzoj1013][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
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 = sqrt( (a1-b1)^2 + (a2-b2)^2 +
… + (an-bn)^2 )
这是一道高斯消元的模板题，我们有n+1个点，设球心的坐标为(x1,x2,x3......,xn)那么每一个点到它的距离都应该是相等的，这样我们可以得到n+1个非线性方程
然后利用n+1个非线性方程，每一个除1以外的非线性方程都可以通过减去第一个非线性方程得到一个线性方程，这样我们就有个n个线性方程
然后把常数项移到等号另一边就可以列出行列式利用高斯消元求解了。
上代码
1 #include<iostream> 2 #include<cstdio> 3 #include<cstdlib> 4 #include<cstring> 5 #define eps 1e-6 6 using namespace std; 7 int n; 8 struct matrix 9 { 10 double d[15]; 11 }s[15],a[15]; 12 double ans[105]; 13 void waps(int x,int y) 14 { 15 matrix tmp; 16 tmp=a[x]; 17 a[x]=a[y]; 18 a[y]=tmp; 19 } 20 double asb(double x) 21 { 22 if(x<0)return -x; 23 return x; 24 } 25 void gauss() 26 { 27 for(int i=1,col=1;i<=n;++i) 28 { 29 int p=i; 30 for(int j=i+1;j<=n;++j) 31 { 32 if(asb(a[j].d[col])-asb(a[p].d[col])>eps) 33 { 34 p=j; 35 } 36 } 37 if(p!=i) 38 waps(p,i); 39 if(a[i].d[col]==0) 40 { 41 col++; 42 i--; 43 continue; 44 } 45 for(int j=i+1;j<=n;++j) 46 { 47 double now=a[i].d[col]/a[j].d[col]; 48 for(int k=col;k<=n+1;++k) 49 { 50 a[j].d[k]*=now,a[j].d[k]-=a[i].d[k]; 51 } 52 } 53 col++; 54 } 55 } 56 void get_ans() 57 { 58 for(int i=n;i>0;--i) 59 { 60 for(int j=n;j>i;--j) 61 { 62 a[i].d[n+1]-=a[i].d[j]*ans[j]; 63 } 64 ans[i]=a[i].d[n+1]/a[i].d[i]; 65 } 66 for(int i=1;i<n;++i) 67 { 68 printf("%.3lf ",ans[i]); 69 } 70 printf("%.3lf",ans[n]); 71 } 72 int main() 73 { 74 scanf("%d",&n); 75 for(int i=0;i<=n;++i) 76 for(int j=1;j<=n;++j) 77 { 78 scanf("%lf",&s[i].d[j]); 79 } 80 for(int i=1;i<=n;++i) 81 for(int j=1;j<=n;++j) 82 { 83 a[i].d[j]=2*s[i].d[j]-2*s[0].d[j]; 84 a[i].d[n+1]+=s[i].d[j]*s[i].d[j]-s[0].d[j]*s[0].d[j]; 85 } 86 gauss(); 87 get_ans(); 88 return 0; 89 }
