洛谷P4035【[JSOI2008]球形空间产生器】题解

传送门

方法：高斯消元

题意简述：给定 \(n\) 个球面上的点的坐标，求球心坐标。

首先，可以肯定的是球面上每个点到球心的距离都是一样的。

则我们看其中一点 \(a_1\) ，它的坐标为 \((a_{1,1},a_{1,2},\cdots,a_{1,n})\) ，设球心为 \(G\) ，坐标为 \((G_1,G_2,\cdots,G_n)\) ，各个点到球心的距离 \(r\) 。

根据题目给我们的公式，可以知道：

\[r=\sqrt{(a_{1,1}-G_1)^2+(a_{1,2}-G_2)^2+\cdots+(a_{1,n}-G_n)^2} \]

\[\therefore r^2=(a_{1,1}-G_1)^2+(a_{1,2}-G_2)^2+\cdots+(a_{1,n}-G_n)^2 \]

将括号拆开，得：

\[r^2=(a_{1,1})^2-2a_{1,1}G_1+(G_1)^2+(a_{1,2})^2-2a_{1,2}G_2+(G_2)^2+\cdots+(a_{1,n})^2-2a_{1,n}G_n+(G_n)^2 \]

移项，得：

\[2(a_{1,1}G_1+a_{1,2}G_2+\cdots+a_{1,n}G_n)+r^2-(G_1)^2-(G_2)^2-\cdots-(G_n)=(a_{1,1})^2+(a_{1,2})^2+\cdots+(a_{1,n})^2 \]

因为 \(a\) 的坐标题目都给我们了， \((a_{1,1})^2+(a_{1,2})^2+\cdots+(a_{1,n})^2\) 的值可以直接求出来，将它的值设为 \(b_1\) 。

对于每一个点 \(a_i\) 我们都有以上的条件，所以 \(r^2-(G_1)^2-(G_2)^2-\cdots-(G_n)\) 可以全部抵消，则可以列出方程组：

\(\begin{cases}a_{1,1}G_1+a_{1,2}G_2+\cdots+a_{1,n}G_n=b_1 \\a_{2,1}G_1+a_{2,2}G_2+\cdots+a_{2,n}G_n=b_2\\\vdots\\a_{n,1}G_1+a_{n,2}G_2+\cdots+a_{n,n}G_n=b_n\\\end{cases}\)

利用高斯消元求解即可。

#include<bits/stdc++.h>
using namespace std;
const double eps=1e-8;
int n;
double zg[102][102],b[102][102];
void gaos(){
	for(int i=1;i<=n;i++){
		int p=i;
		for(int j=i+1;j<=n;j++){
			p=(abs(zg[j][i])>abs(zg[p][i])?j:p);
		}
		if(abs(zg[p][i])<eps){
			return ;
		}
		if(p!=i){
			for(int j=i;j<=n+1;j++){
				swap(zg[i][j],zg[p][j]);
			}
		}
		for(int j=i+1;j<=n;j++){
			if(i==j) continue;
			double x=zg[j][i]/zg[i][i];
			for(int k=i;k<=n+1;k++){
				zg[j][k]=zg[j][k]-zg[i][k]*x;
			}
		}
	}
	for(int i=n;i;i--){
		for(int j=i+1;j<=n;j++){
			zg[i][n+1]=zg[i][n+1]-zg[j][n+1]*zg[i][j];
		}
		zg[i][n+1]=zg[i][n+1]/zg[i][i];
	}
}
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++){
			zg[i][j]=2*(b[i][j]-b[i+1][j]);
			zg[i][n+1]+=b[i][j]*b[i][j]-b[i+1][j]*b[i+1][j];
		}
	}
	gaos();
	for(int i=1;i<=n;i++){
		printf("%.3lf ",zg[i][n+1]);
	}
	return 0;
}