洛谷4035 JSOI2008球形空间产生器 (列柿子+高斯消元)

题目链接


qwq
首先看到这个题,感觉就应该从列方程入手。

我们设给定的点的坐标矩阵是\(x\),然后球心坐标\(a_1,a_2....a_n\)

根据欧几里得距离公式,对于一个\(n维空间\)的第\(i\)个点,他距离球心的距离可以表示为\[\sum_{j=1}^n (x_{ij}-a[j])^2 = r^2 \]

通过\(n+1\)个点,我们可以轻松列出来\(n+1\)个方程,但是可惜不是线性。我们考虑该怎么优化这个过程

考虑到相邻的两个方程的右边都是相等的,我们不妨将相邻两个方程进行减法,得到n个新的方程,假设\(i和i+1\)进行减法
\[\sum_{j=1}^n (x_{ij}-a[j])^2 - \sum_{j=1}^n (x_{i+1j}-a[j])^2 = 0\]
\[\sum_{j=1}^n x_{ij}^2 - \sum_{j=1}^n x_{i+1j}^2 + 2\times \sum(x_{i+1j}-x_{ij})*a_j = 0\]

其中\(x\)是已知量,我们稍微移项,就可以得到最后的柿子啦
\[\sum(x_{i+1j}-x_{ij})*a_j = \frac{\sum_{j=1}^n x_{i+1j}^2 - \sum_{j=1}^n x_{ij}^2}{2}\]

那么经过一番变化,现在我们有n个未知数,n个方程,正好可以解出来每个未知数的解。

直接上高斯消元即可

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
#include<queue>
#include<map>
#include<set>
#define mk make_pair
#define ll long long
using namespace std;
inline int read()
{
  int x=0,f=1;char ch=getchar();
  while (!isdigit(ch)) {if (ch=='-') f=-1;ch=getchar();}
  while (isdigit(ch)) {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
  return x*f;
}
const int maxn = 510;
const double eps = 1e-8;
double a[maxn][maxn];
int n,m;
double x[maxn];
double b[maxn][maxn];
double sum[maxn];
void gauss()
{
    int k=1;
    for (int i=1;i<=n;i++)
    {
        int now = k;
        while (now<=n && fabs(a[now][i])<=eps) now++;
        if (now==n+1) continue;
        for (int j=1;j<=n+1;j++) swap(a[now][j],a[k][j]);
        for (int j=1;j<=n;j++)  
        {
             if (j!=k)
             {
                 double t = a[j][i]/a[k][i];
                 for(int p=1;p<=n+1;p++) a[j][p]=a[j][p]-a[k][p]*t; 
             }
        }
        ++k;
    }
    for (int i=n;i>=1;i--)
    {
        for (int j=i+1;j<=n;j++) a[i][n+1]-=a[i][j]*x[j];
        x[i]=a[i][n+1]/a[i][i];
    }
} 
int main()
{
  n=read();
  for (int i=1;i<=n+1;i++)
  {
     for (int j=1;j<=n;j++)
       scanf("%lf",&b[i][j]),sum[i]+=b[i][j]*b[i][j];
  } 
  for (int i=1;i<=n;i++)
  {
     for (int j=1;j<=n;j++)
       a[i][j]=b[i+1][j]-b[i][j];
     a[i][n+1]=(sum[i+1]-sum[i])/2;
  }
  gauss();
  for (int i=1;i<=n;i++) printf("%.3lf ",x[i]);
  return 0;
}

转载于:https://www.cnblogs.com/yimmortal/p/10151318.html

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值