n维球上的任意一点到球心距离相等，故设球心坐标为(x1,x2,...,xn)
则有公式∑(ai,j-xj)2=C
设法消去平方项x2，考虑相邻两项相减
可将公式化为Σ2(ai,j-ai+1,j)xj=Σ(ai,j2-ai+1.j2)
此时就可以构造线性方程组了
#include<iostream> #include<cstdio> #include<algorithm> #include<cstring> #include<cmath> #include<queue> using namespace std; double eps=1e-8; double a[1010][1010],b[1010][1010]; int n; int flag=0; bool pd(int k,int x,int y){ if (fabs(fabs(a[x][k])-fabs(a[y][k]))>eps) return fabs(a[x][k])>fabs(a[y][k]); for (int i=k+1;i<=n;i++){ if (fabs(fabs(a[x][i])-fabs(a[y][i]))>eps) return fabs(a[x][i])<fabs(a[y][i]); } return 0; } void print() { int i,j; for(i=1;i<=n;i++) { for(j=1;j<=n+1;j++) { printf("%f ,",a[i][j]); } printf("\n"); } } void swaprow(int now,int maxrow) { double tmp; for(int i=now;i<=n+1;i++) { tmp=a[now][i]; a[now][i]=a[maxrow][i]; a[maxrow][i]=tmp; } return; } void Gauss() { double tmp; int maxrow; for(int j=1;j<=n;j++) { maxrow=j; for(int i=j;i<=n;i++) { if(pd(j,i,maxrow))maxrow=i; maxrow=i; }//选择主元 if(maxrow!=j) { swaprow(j,maxrow); } for(int i=j+1;i<=n;i++) { tmp=a[i][j]/a[j][j]; for(int k=j;k<=n+1;k++) { a[i][k]-=a[j][k]*tmp; } } //print(); //cout<<endl; } for(int i=n;i>=1;i--) { for(int j=i+1;j<=n;j++) a[i][n+1]-=a[i][j]*a[j][n+1]; if(a[i][n+1]==0&&a[i][i]==0){ flag=1; } if(a[i][n+1]!=0&&a[i][i]==0){ flag=2; return; } a[i][n+1]/=a[i][i]; } return ; } int main() { cin>>n; for(int i=1;i<=n+1;i++) { for(int j=1;j<=n;j++) { cin>>b[i][j]; } } for(int i=1;i<=n;i++) { for(int j=1;j<=n;j++) { a[i][j]=2*(b[i][j]-b[i+1][j]); a[i][n+1]+=(b[i][j]*b[i][j]-b[i+1][j]*b[i+1][j]); } } Gauss(); for(int i=1;i<=n;i++) { if(fabs(a[i][n+1])<eps)a[i][n+1]=0; printf("%.3lf ",a[i][n+1]); } return 0; }