【洛谷P4035】【BZOJ1013】球形空间产生器【高斯消元】

题目大意:

题目链接:
洛谷:https://www.luogu.org/problem/P4035
BZOJ:https://www.lydsy.com/JudgeOnline/problem.php?id=1013
给出nn维中的n+1n+1个在nn维球体球面上的点,求出这个球的球心坐标。


思路:

题目下方有如下说明:
距离:设两个n为空间上的点A,BA, B的坐标为(a1,a2,...,an),(b1,b2,...bn)(a_1,a_2,...,a_n),(b_1,b_2,...b_n),则AB的距离定义为:dist=(a1b1)2+(a2b2)2+...+(anbn)2dist=\sqrt{(a_1-b_1)^2+(a_2-b_2)^2+...+(a_n-b_n)^2}
那么我们就是要在nn维平面上找一个点(x1,x2,...,xn)(x_1,x_2,...,x_n),满足
j=1n(a1,jxi)2=j=1n(a2,jxi)2=...=j=1n(an+1,jxi)2\sum^{n}_{j=1}(a_{1,j}-x_i)^2=\sum^{n}_{j=1}(a_{2,j}-x_i)^2=...=\sum^{n}_{j=1}(a_{n+1,j}-x_i)^2
把等号左右的式子相减,以第一个等号维例
j=1n(a1,j2+xi22a1,jxi(a2,j2+xi22a2,jxi))\sum^{n}_{j=1}(a_{1,j}^2+x_i^2-2a_{1,j}x_i-(a_{2,j}^2+x_i^2-2a_{2,j}x_i))
j=1n(a1,j2a2,j22xi(a1,ja2j))=0\sum^{n}_{j=1}(a_{1,j}^2-a_{2,j}^2-2x_i(a_{1,j}-a_{2_j}))=0
j=1n2(a1,ja2,j)xj=j=1n(a1,j2a2,j2)\sum^{n}_{j=1}2(a_{1,j}-a_{2,j})x_j=\sum^{n}_{j=1}(a_{1,j}^2-a_{2,j}^2)
然后就可以构造出增广矩阵
[2(a1,1a2,1)2(a1,2a2,2)2(a1,na2,n)j=1n(a1,j2a2,j2)2(a2,1a3,1)2(a2,2a3,2)2(a2,na3,n)j=1n(a2,j2a3,j2)2(an,1an+1,1)2(an,2an+1,2)2(an+1,nan,n)j=1n(an,j2an+1,j2)]\begin{bmatrix} &2(a_{1,1}-a_{2,1}) &2(a_{1,2}-a_{2,2}) &\cdots &2(a_{1,n}-a_{2,n}) &| & \sum^{n}_{j=1}(a_{1,j}^{2}-a_{2,j}^2)\\ &2(a_{2,1}-a_{3,1}) &2(a_{2,2}-a_{3,2}) &\cdots &2(a_{2,n}-a_{3,n}) &| & \sum^{n}_{j=1}(a_{2,j}^{2}-a_{3,j}^2)\\ &\vdots & \vdots & \ddots & \vdots & | & \vdots\\ &2(a_{n,1}-a_{n+1,1}) &2(a_{n,2}-a_{n+1,2}) & \cdots & 2(a_{n+1,n}-a_{n,n}) & | &\sum^{n}_{j=1}(a_{n,j}^{2}-a_{n+1,j}^2) \end{bmatrix}
md打的累死我了
然后套高斯消元的模板就好了


代码:

#include <cmath>
#include <cstdio>
#include <algorithm>
using namespace std;

const int N=15;
double a[N][N],b[N],c[N][N];
int n;

int main()
{
	scanf("%d",&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++)
		for (int j=1;j<=n;j++)
		{
			b[i]+=a[i][j]*a[i][j]-a[i+1][j]*a[i+1][j];
			c[i][j]=(a[i][j]-a[i+1][j])*2.0;
		}
	for (int i=1;i<=n;i++)
	{
		for (int j=i;j<=n;j++)
			if (fabs(c[j][i])>1e-8)
			{
				for (int k=1;k<=n;k++) swap(c[j][k],c[i][k]);
				swap(b[j],b[i]);
				break;
			}
		for (int j=1;j<=n;j++)
			if (i!=j)
			{
				double rate=c[j][i]/c[i][i];
				for (int k=1;k<=n;k++) c[j][k]-=c[i][k]*rate;
				b[j]-=b[i]*rate;
			}
	}
	for (int i=1;i<=n;i++)
		printf("%0.3lf ",b[i]/c[i][i]);
	return 0;
}
posted @ 2019-08-16 09:57  全OI最菜  阅读(121)  评论(0编辑  收藏  举报