BZOJ4689 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.385379805

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

乍一看以为是差值

~~细细一想没有给出x？~~

定睛观察发现给出了点值和对应的x

然后笑嘻嘻的开始写差值

结果又细细一想

woc，这不是个高斯消元吗？

关于系数的高斯消元！！！

然后就裸了

就是被精度卡了一下

还有一个特别坑的是题目中的是从x=0开始给的

#include<bits/stdc++.h>

using namespace std;

#define N 10

#define eps 1e-3

double p[N],a[N][N];

void Guass(int n,double g[N][N]){

    for(int i=1;i<=n;i++){

        int p=0;

        for(int j=i;j<=n;j++)

            if(g[j][i]>eps){p=j;break;}

        for(int j=1;j<=n+1;j++)swap(g[p][j],g[i][j]);

        double d=g[i][i];

        for(int j=1;j<=n+1;j++)g[i][j]/=d;

        for(int j=1;j<=n;j++){

            if(i==j)continue;

            double w=g[j][i];

            for(int k=1;k<=n+1;k++)

                g[j][k]-=g[i][k]*w;

        }

    }

}

double fast_pow(double a,int b){

    double ans=1.0;

    while(b){

        if(b&1)ans*=a;

        b>>=1;

        a*=a;

    }

    return ans;

}

bool check(int pos,int n){

    static bool vis[N];

    memset(vis,0,sizeof(vis));

    vis[pos]=1;

    int cnt=0;

    for(int i=0;i<=n+1;i++){

        if(i==pos)continue;

        cnt++;vis[i]=1;

        for(int j=1;j<=n;j++)a[cnt][j]=fast_pow((double)i,j-1);

        a[cnt][n+1]=p[i];

        if(cnt==n)break;

    }

    Guass(n,a);

    for(int i=0;i<=n+1;i++){

        if(vis[i])continue;

        double res=0;

        for(int j=1;j<=n;j++)res+=a[j][n+1]*fast_pow((double)i,j-1);

        if(fabs(res-p[i])>eps)return 0;

    }

    return 1;

}

int main(){

    int n;

    while(1){

        scanf("%d",&n);

        if(!n)return 0;

        n++;//实际上是n+1个变量

        for(int i=0;i<=n+1;i++)scanf("%lf",&p[i]);//f[0]~f[n+1]

        int ans=0;

        for(int i=0;i<=n+1;i++)

            if(check(i,n)){ans=i;break;}

        printf("%d\n",ans);

    }

    return 0;

}

巴特西

BZOJ4689 Find the Outlier 【高斯消元】*