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

2019-05-14 16:50:31


本题的难点在于建方程。

设球半径为 $r$ ,球心为 $a1$ , $a2$ ,......, $an$ ,坐标分别为 $P_{1,1}$ , $P_{1,2}$ ,......, $P_{n+1,n}$

考虑 $P_1$ ,有方程 $\sum_{k=1}^{k<=n}(P_{1,k}-a_k)^2=r^2$

展开得 $\sum_{k=1}^{k<=n}(P_{1,k})^2-2(P_{1,k})(a_k)+(a_k)^2=r^2$

同理,考虑 $P_2$ ,有方程 $\sum_{k=1}^{k<=n}(P_{2,k}-a_k)^2=r^2$

展开得 $\sum_{k=1}^{k<=n}(P_{2,k})^2-2(P_{2,k})(a_k)+(a_k)^2=r^2$

两个展开式相减,得 $\sum_{k=1}^{k<=n}[(P_{1,k})^2-2(P_{1,k})(a_k)+(a_k)^2-(P_{2,k})^2+2(P_{2,k})(a_k)-(a_k)^2]=0$

即 $\sum_{k=1}^{k<=n}[(P_{1,k})^2-2(P_{1,k})(a_k)-(P_{2,k})^2+2(P_{2,k})(a_k)]=0$

$\sum_{k=1}^{k<=n}[-2(P_{1,k})(a_k)+2(P_{2,k})(a_k)]=\sum_{k=1}^{k<=n}[-(P_{1,k})^2+(P_{2,k})^2]$

此时原方程次数为1,可以求解。

上代码:

#include <bits/stdc++.h>
using namespace std;
int n;
double a[15][15];
long double mat[15][15],ans[15];
inline double sqr(double x) {return x * x;}
const long double EPS = 1e-10;
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++) 
        {
            mat[i][j] = 2 * (a[i + 1][j] - a[i][j]);//未知数系数
            mat[i][n + 1] += sqr(a[i + 1][j]) - sqr(a[i][j]);//常数(
        }
    for (int i = 1;i <= n;i++)//高斯消元
    {
        auto mx = abs(mat[i][i]);int pos = i;
        for (int j = i + 1;j <= n;j++)//寻找系数绝对值最大的列提高精度
            if (abs(mat[j][i]) > mx)
            {
                mx = abs(mat[j][i]);pos = j;
            }
        swap(mat[i],mat[pos]);
        for (int j = i + 1;j <= n;j++)//消元
        {
            auto d = mat[j][i] / mat[i][i];
            for (int k = i + 1;k <= n + 1;k++) mat[j][k] -= mat[i][k] * d;
            mat[j][i] = 0;
        }
    }
    for (int i = n;i > 0;--i)//回代
    {
        for (int j = i + 1;j <= n;j++)
            mat[i][n + 1] -= ans[j] * mat[i][j];
        ans[i] = mat[i][n + 1] / mat[i][i];
    }
    for (int i = 1;i <= n;i++) printf("%.3Lf ",ans[i]);
    return 3;
}