题目链接:https://www.acwing.com/problem/content/209/
给出n维空间中的n维圆上的n+1个点,要求求这个圆的圆心,保证一定有解,这样一共有n+1个n圆二次方程,差分之后变成n个n元方程,保证一定有解那么系数矩阵一定是满秩的,通过高斯消元可以求解
代码:
#include<iostream> #include<cstdio> #include<algorithm> #include<math.h> using namespace std; #define maxn 20 double a[maxn][maxn],b[maxn],c[maxn][maxn]; int n; int main(){ cin>>n; for(int i=1;i<=n+1;i++) for(int j=1;j<=n;j++) scanf("%lf",&a[i][j]); for(int i=1;i<=n;i++){ //构造n个差分方程, for(int j=1;j<=n;j++){ c[i][j]=2*(a[i][j]-a[i+1][j]); b[i]+=a[i][j]*a[i][j]-a[i+1][j]*a[i+1][j]; } } //高斯消元 for(int i=1;i<=n;i++){ for(int j=i;j<=n;j++){//只交换比i大的方程j if(fabs(c[j][i])>1e-8){//系数不为0 for(int k=1;k<=n;k++)swap(c[j][k],c[i][k]);//交换两个方程 swap(b[i],b[j]); } } for(int j=1;j<=n;j++){//第j行减去第i行的rate倍 if(j==i)continue;//将其他方程的对应项全部变成0 double rate=c[j][i]/c[i][i]; for(int k=i;k<=n;k++){//i之前的全都是0,不用处理 c[j][k]-=c[i][k]*rate; } b[j]-=b[i]*rate; } } for(int i=1;i<=n;i++)b[i]/=c[i][i],c[i][i]=1;//将x_i的系数变为1, for(int i=1;i<n;i++)printf("%.3lf ",b[i]); printf("%.3lf ",b[n]); }