[BZOJ]4689: Find the Outlier 高斯消元

Description

Abacus教授刚刚完成了一个制作数表的计算引擎的设计。它被设计用于同时计算一个多项式在许多点的取值。例如对于多项式 f(x)=x^2+2x+1 ,一种可能的计算结果是 f(0)=1,f(1)=4,f(2)=9.f(3)=16,f(4)=25 。不幸的是,引擎存在一个故障使得计算出的值总有一个是错的,例如对于上述多项式,它可能输出 1,4,12,16,25 而不是 1,4,9,16,25 。请你帮教授找出发生故障的是哪个点值。
Input

输入包含多组测试数据。
每组数据第一行包含一个正整数 d 表示多项式的度数,即多项式最高次项的项数,保证 d≤5 。
接下来 d+3 行,每行一个实数,第 i 行表示输出的 f(i) 的值,保证-100.0≤f(i)≤100.0 。
你可以认为恰好只有一个点值出故障,且与实际值的误差超过 1.0 。
由于不可避免的误差,其他数字与精确值的误差不超过10^(-6) 。
输入以一个零作为结束。
Output

对于每组数据,输出一个非负整数 i 表示 f(i) 的值发生故障。
Sample Input

2

1.0

4.0

12.0

16.0

25.0

1

-30.5893962764

5.76397083962

39.3853798058

74.3727663177

4

42.4715310246

79.5420238202

28.0282396675

-30.3627807522

-49.8363481393

-25.5101480106

7.58575761381

5

-21.9161699038

-48.469304271

-24.3188578417

-2.35085940324

-9.70239202086

-47.2709510623

-93.5066246072

-82.5073836498

0
Sample Output

2

1

1

6

题解:高斯消元

这道水题不知道为什么这么少人做……因为数据范围极小,直接暴力枚举哪个是错的,然后用高斯消元check即可。注意eps要开大,不然会WA。

代码:

#include<cstdio>
#include<cstring>
#include<cmath>
#include<iostream>
#include<algorithm>
using namespace std;
const double eps=1e-3;
double f[100],s[100][100],ans[100];
double p(double x,int y)
{
    double re=1.0;
    for(int i=1;i<=y;i++)re*=x;
    return re;
}
int n;
bool check(int x)
{
    bool mark[10];
    memset(mark,false,sizeof(mark));
    int cnt=0;mark[x]=true;
    for(int i=1;i<n+2;i++)
    {
        if(!mark[i])
        {
            mark[i]=true;
            cnt++;
            for(int j=1;j<=n;j++)
            s[cnt][j]=p((double)i,j-1);
            s[cnt][n+1]=f[i];
        }
        if(cnt==n)break;
    }
    for(int i=1;i<n;i++)
    {
        double t1=s[i][i];
        for(int j=1;j<=n+1;j++) s[i][j]/=t1;
        for(int j=i+1;j<=n;j++) 
        {
            double t2=s[j][i];
            for(int k=1;k<=n+1;k++) s[j][k]-=(s[i][k]*t2);
        }
    }
    ans[n]=s[n][n+1]/s[n][n];
    for(int i=n-1;i>0;i--)
    {
        double num=s[i][n+1];
        for(int j=n;j>i;j--) num-=(s[i][j]*ans[j]);
        ans[i]=num/s[i][i];
    }
    for(int i=0;i<n+2;i++)
    if(!mark[i])
    {
        double t=0.0;
        for(int j=0;j<n;j++)
        t+=ans[j+1]*p((double)i,j);
        if(abs(t-f[i])>=eps)return false;
    }
    return true;
}
int main()
{
    while(1)
    {
        scanf("%d",&n);
        if(!n)break;
        for(int i=0;i<n+3;i++)scanf("%lf",&f[i]);
        n++;
        bool tf=false;
        for(int i=1;i<n+2;i++)
        if(check(i)){printf("%d\n",i);tf=true;break;}
        if(!tf)puts("0");
    }
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值