AcWing 207. 球形空间产生器

传送门

思路：

        设球心坐标为（x1，x2，...，xn），有 ，由此我们可以列出N+1个二次方程，我们可以对前后两个方程做差，来得到N个一次方程，同时可以消掉常数C，第i个方程即

那么我们就可以直接采用高斯消元，解出圆心的坐标。

代码： 

#include<bits/stdc++.h>
#include<unordered_map>
using namespace std;
typedef long long LL;
typedef unsigned long long ULL;
typedef pair<int, int> PII;
typedef pair<double, double> PDD;
//#define int LL
#define inf 0x3f3f3f3f
#define INF 0x3f3f3f3f3f3f3f3f
#define IOS ios::sync_with_stdio(0),cin.tie(0),cout.tie(0)
#pragma warning(disable :4996)
const int maxn = 2000100;
const double eps = 1e-8;
const LL MOD = 998244353;

double a[20][20], b[20], c[20][20];//b与c构成增广矩阵
int N;
double A[20][20];

int gauss()
{
	int ans = N;//主元个数
	for (int i = 1; i <= N; i++)
	{
		int temp = i;
		for (int j = i; j <= N; j++)
		{
			if (fabs(a[j][i]) > fabs(a[temp][i]))
				temp = j;
		}
		if (fabs(a[temp][i]) < eps)//当前列无主元
		{
			ans--;
			continue;
		}
		for (int k = 1; k <= N; k++)//第i列系数不为0的行换到第i行
			swap(a[i][k], a[temp][k]);
		swap(b[i], b[temp]);
		double div1 = a[i][i];
		for (int j = i; j <= N; j++)
			a[i][j] /= div1;
		b[i] /= div1;
		for (int j = 1; j <= N; j++)
		{
			if (i != j)
			{
				double div2 = a[j][i] / 1.0;
				for (int k = i; k <= N; k++)
					a[j][k] -= a[i][k] * div2;
				b[j] -= b[i] * div2;
			}
		}
	}

	return ans;
}

void solve()
{
	for (int i = 1; i <= N; i++)
	{
		b[i] = 0;
		for (int j = 1; j <= N; 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 = 1; j <= N; j++)
			a[i][j] = 2.0 * (A[i][j] - A[i + 1][j]);
	}
	gauss();
	for (int i = 1; i <= N; i++)
		printf("%.3lf ", b[i]);
	putchar('\n');
}

int main()
{
	IOS;
	scanf("%d", &N);
	for (int i = 1; i <= N + 1; i++)
	{
		for (int j = 1; j <= N; j++)
			scanf("%lf", &A[i][j]);
	}
	solve();

	return 0;
}