令球心坐标为x1,x2...xn,假设当前第i个点坐标为a1,a2...,an,第i+1个点坐标为b1,b2...,bn,则由半径相等可得:
(a1-x1)^2+(a2-x2)^2+...+(an-xn)^2=(b1-x1)^2+(b2-x2)^2+...+(bn-xn)^2
化简可得:
2(a1-b1)x1+2(a2-b2)x2+...+2(an-bn)xn=(a1^2+a2^2+...+an^2-b1^2-b2^2-...-b3^2)
如此可得到n个一元n次方程组,用最简单的高斯消元搞一搞就好了。
By:大奕哥
1 #include<bits/stdc++.h> 2 using namespace std; 3 const int N=15; 4 double a[N][N],x[N][N],sum[N],ans[N]; 5 int n; 6 void Gauss() 7 { 8 for(int i=1;i<=n;++i) 9 { 10 int t=i; 11 for(int j=i+1;j<=n;++j)if(a[j][i]>a[t][i])t=j; 12 if(t!=i)for(int j=1;j<=n+1;j++)swap(a[i][j],a[t][j]); 13 for(int j=i+1;j<=n;++j) 14 { 15 double tmp=a[j][i]/a[i][i]; 16 for(int k=i;k<=n+1;++k)a[j][k]-=tmp*a[i][k]; 17 } 18 } 19 for(int i=n;i;--i) 20 { 21 ans[i]=a[i][n+1]/a[i][i]; 22 for(int j=i-1;j;--j) 23 { 24 a[j][n+1]-=a[j][i]*ans[i]; 25 } 26 } 27 } 28 int main() 29 { 30 scanf("%d",&n); 31 for(int i=1;i<=n+1;++i) 32 { 33 for(int j=1;j<=n;++j) 34 scanf("%lf",&x[i][j]),sum[i]=sum[i]+x[i][j]*x[i][j]; 35 } 36 for(int i=1;i<=n;++i) 37 { 38 for(int j=1;j<=n;++j)a[i][j]=2*(x[i][j]-x[i+1][j]); 39 a[i][n+1]=sum[i]-sum[i+1]; 40 } 41 Gauss(); 42 for(int i=1;i<n;++i)printf("%.3lf ",ans[i]); 43 printf("%.3lf\n",ans[n]); 44 return 0; 45 }