68 lines
1.7 KiB
C++
68 lines
1.7 KiB
C++
|
#include<stdio.h>
|
||
|
#include<math.h>
|
||
|
// using newton interpolation
|
||
|
int main(){
|
||
|
int n,m;
|
||
|
double x[100],y[100];
|
||
|
double t[100];
|
||
|
double xt;
|
||
|
scanf("%d",&n);
|
||
|
scanf("%d",&m);
|
||
|
for(int i=0;i<n;i++){
|
||
|
scanf("%lf %lf",&x[i],&y[i]);
|
||
|
for(int j=0;j<i;j++){
|
||
|
// if points are too close, just keep one of them.
|
||
|
if(fabs(x[i]-x[j])<1e-6){
|
||
|
i--;
|
||
|
n--;
|
||
|
continue;
|
||
|
}
|
||
|
}
|
||
|
}
|
||
|
// newton interpolation
|
||
|
for(int i=0;i<n;i++){
|
||
|
for(int j=i-1;j>=0;j--){
|
||
|
double p=1;
|
||
|
for(int k=j-1;k>=0;k--){
|
||
|
p=p*(x[i]-x[k]);
|
||
|
}
|
||
|
y[i]=y[i]-p*t[j];
|
||
|
}
|
||
|
double p=1;
|
||
|
for(int k=i-1;k>=0;k--){
|
||
|
p=p*(x[i]-x[k]);
|
||
|
}
|
||
|
/*
|
||
|
if value is precise enough then higher terms are not necessary (sometimes even
|
||
|
troublesome, when x[k] is too close to x[0], x[1],...,x[k-1] (distance < 1)
|
||
|
then p=(x[k]-x[0])*(x[k]-x[1])*...*(x[k]-x[k-1]) is too small, this will result
|
||
|
in a super large t[k] even when the residual error of a k-1 degree polynomial
|
||
|
is small enough.)
|
||
|
*/
|
||
|
if(fabs(y[i])<1e-9){
|
||
|
y[i]=0;
|
||
|
}
|
||
|
t[i]=y[i]/p;
|
||
|
}
|
||
|
int k=0;
|
||
|
for(int i=0;i<n;i++){
|
||
|
// if a coefficient is too small, then regard it as zero.
|
||
|
if(fabs(t[i])>1e-6){
|
||
|
k=i;
|
||
|
}
|
||
|
}
|
||
|
// degree of the polynomial
|
||
|
printf("%d\n",k);
|
||
|
n=k+1;
|
||
|
for(int i=0;i<m;i++){
|
||
|
scanf("%lf",&xt);
|
||
|
double yt=0;
|
||
|
double p=1;
|
||
|
for(int j=0;j<n;j++){
|
||
|
yt+=t[j]*p;
|
||
|
p*=(xt-x[j]);
|
||
|
}
|
||
|
printf("%lf\n",yt);
|
||
|
}
|
||
|
return 0;
|
||
|
}
|