牛顿插值法的C语言编程

合集下载

Newton插值的C++实现

Newton插值的C++实现

Newton插值的C++实现Newton(⽜顿)插值法具有递推性,这决定其性能要好于Lagrange(拉格朗⽇)插值法。

其重点在于差商(Divided Difference)表的求解。

步骤1. 求解差商表,这⾥采⽤⾮递归法(看着挺复杂挺乱,这⾥就要⾃⼰动笔推⼀推了,闲了补上其思路),这样,其返回的数组(指针)就是差商表了,/** 根据插值节点及其函数值获得差商表* 根据公式⾮递归地求解差商表* x: 插值节点数组* y: 插值节点处的函数值数组* lenX: 插值节点的个数* return: double类型数组*/double * getDividedDifferenceTable(double x[], double y[], int lenX) {double *result = new double[lenX*(lenX - 1) / 2];for (int i = 0; i < lenX - 1; i++) {result[i] = (y[i + 1] - y[i]) / (x[i + 1] - x[i]);}int step = lenX - 1; // 增加步长int yindex = lenX - 1; // 分⼦的基准值,线性减速递增int xgap = 2; // 分母的间距,每次+1int xstep = 2;while (step >= 1) {for (int i = 0; i < step - 1; i++) {result[yindex + i] = (result[yindex - step + i + 1] - result[yindex - step + i]) / (x[xstep + i] - x[xstep + i - xgap]);}yindex += (step - 1);xstep++;step--;xgap++;}return result;}步骤 2. 取得差商表每⼀列的第⼀个作为基函数的系数/*** 从差商表中获取⼀定阶数的某个差商* dd: 差商表* len: 差商表的长度* rank: 所取差商的阶数* index: rank阶差商的第index个差商* return: 返回double类型所求差商*/double getDividedDifference(double dd[], int len, int rank, int index) {// 根据差商表的长度求解差商的最⾼阶数// 由于差商表是三⾓形的,因此有规律可循,解⽅程即可int rankNum = (int)(0.5 + sqrt(2 * len + 0.25) - 1);printf("%d 阶差商 \n", rank);int pos = 0;int r = 1;// 根据n+(n+1)+...+(n-rank)求得rank阶差商在差商表数组中的索引while (rank > 1){// printf("geting dd's index, waiting...\n");pos += rankNum;rank--;rankNum--;}pos += index; // 然后加上偏移量,意为rank阶差商的第index个差商return dd[pos];}步骤3. Newton 插值主流程/** Newton Interpolating ⽜顿插值主要代码* x: 插值节点数组* y: 插值节点处的函数值数组* lenX: 插值节点个数(数组x的长度)* newX: 待求节点的数组* lennx: 待求节点的个数(数组lenX的长度)*/double *newtonInterpolation(double x[], double y[], int lenX, double newX[], int lennx) {// 计算差商表double *dividedDifferences = getDividedDifferenceTable(x, y, lenX);//double *result = new double[lennx];// 求差商表长度int ddlength = 0;for (int i = 1; i < lenX; i++){ddlength += i;}// 打印差商表信息printf("=======================差商表的信息=======================\n");printf("差商个数有:%d, 待插值节点个数:%d\n =======================差商表=======================", ddlength, lennx); int ifnextrow = 0; // 控制打印换⾏for (int i = 0; i < ddlength; i++) {if (ifnextrow == (i)*(i + 1) / 2)printf("\n");printf("%lf, ", dividedDifferences[i]);ifnextrow++;}printf("\n");// 计算Newton Interpolating 插值double ddi = 0;double coef = 1;for (int i = 0; i < lennx; i += 2) {coef = (double)1;result[i] = y[i];printf("=============打印计算过程=============\n");for (int j = 1; j < lenX -1; j++) {coef *= (newX[i] - x[j - 1]);ddi = getDividedDifference(dividedDifferences, ddlength, j, 0);printf("取得差商: %lf\n", ddi);result[i] = result[i] + coef * ddi;printf("计算result[%d] + %lf * %lf", i, coef, ddi);printf("获得结果result %d: %lf\n", i, result[i]);}printf("result[%d] = %lf\n", i, result[i]);// =======================选做题:求误差==========================printf("求解截断误差所需差商:%lf\n", getDividedDifference(dividedDifferences, ddlength, lenX-1, 0));// printf("求解截断误差所需多项式:%lf\n", coef);printf("求解截断误差所需多项式的最后⼀项:%lf - %lf\n", newX[i] , x[lenX - 2]);result[i + 1] = getDividedDifference(dividedDifferences, ddlength, lenX - 1, 0)*coef*(newX[i] - x[lenX - 2]);// ===============================================================}if (dividedDifferences) {delete[] dividedDifferences;}return result;}步骤4. 主函数⽰例int main(){std::cout << "Hello World!\n";printf("Newton插值!\n");double x[] = {0.40, 0.55, 0.65, 0.80, 0.90, 1.05};//double x[] = { 1.5, 1.6, 1.7 };int lenX = sizeof(x) / sizeof(x[0]);double y[] = {0.41075, 0.57815, 0.69675, 0.88811, 1.02652, 1.25382};//double y[] = { 0.99749, 0.99957, 0.99166 };double newx[] = {0.596};//double newx[] = { 1.609 };// 计算⽜顿插值结果和截断误差double *result = newtonInterpolation(x, y, lenX, newx, 1);printf("=====================最终结果=====================\n");printf("求得sin(1.609)的近似值为:%lf\n", *result);printf("截断误差:%.10lf\n", *(result + 1));if (result) {delete[] result;}return0;}。

拉格朗日插值牛顿插值C语言实验报告

拉格朗日插值牛顿插值C语言实验报告

实验报告:数学与统计学系信息与计算科学专业实验报告一、题目1、上机作业题程序12、上机作业题程序2二、算法1、Lagrange 插值//输入被插值点的数目POINT;int main(){int n;inti,j;POINT points[MAX_N+1];double diff[MAX_N+1];doublex,tmp=0,lagrange=0,tx,ty;printf("\nInput n value:");scanf("%d",&n);if(n>MAX_N){printf("The input n is larger thenMAX_N,please redefine the MAX_N.\n");return 1;}if(n<=0){printf("Please input a number between 1 and %d\n",MAX_N);return 1;}//输入被插值点printf("Now input the (x_i,y_i),i=0,...,%d:\n",n);for(i=0;i<=n;i++)scanf("%lf%lf",&points[i].x,&points[i].y);printf("Now input the x value:"); //输入计算Lagrange插值多项式的x值scanf("%lf",&x);for(i=0;i<=n;i++){diff[i]=0;tx=1;ty=1;for(j=0;j<=n;j++){if(i!=j){tx=tx*(x-points[j].x);ty=ty*(points[i].x-points[j].x);}}diff[i]=tx/ty;}for(i=0;i<=n;i++){tmp=points[i].y*diff[i];printf("%f",tmp);lagrange+=tmp;}printf("lagrange(%f)=%f\n",x,lagrange);return 0;}2、Newton 插值//输入被插值点的数目POINT;int main(){ int n;inti,j;POINT points[MAX_N+1];double diff[MAX_N+1];doublex,tmp,newton=0;printf("\nInput n value: ");scanf("%d",&n);if (n>MAX_N){printf("The input n is larger thenMAX_N,please redefine the MAX_N.\n");return 1;}if(n<=0){printf("Please input a number between 1 and %d\n",MAX_N);// getch(); return 1;}//输入被插值点printf("Now input the (x_i,y_i),i=0,...,%d:\n",n);for (i=0;i<=n;i++)scanf("%lf%lf",&points[i].x,&points[i].y);printf("Now input the x value: ");//输入计算Newton插值多项式的x值scanf("%lf",&x);for (i=0;i<=n;i++)diff[i]=points[i].y;for (i=0;i<n;i++){for (j=n;j>i;j--){diff[j]=(diff[j]-diff[j-1])/(points[j].x-points[j-1-i].x);}//计算f(x_0,…,x_n)的差商}tmp=1;newton=diff[0];for(i=0;i<n;i++){tmp=tmp*(x-points[i].x);newton=newton+tmp*diff[i+1];}printf("newton(%f)=%f\n",x,newton);return 0;}三、C程序1、Lagrange 插值#include <stdio.h>#define MAX_N 20typedefstructtagPOINT{double x;double y;}POINT;int main(){int n;inti,j;POINT points[MAX_N+1];double diff[MAX_N+1];doublex,tmp=0,lagrange=0,tx,ty;printf("\nInput n value:");scanf("%d",&n);if(n>MAX_N){printf("The input n is larger thenMAX_N,please redefine the MAX_N.\n");return 1;}if(n<=0){printf("Please input a number between 1 and %d\n",MAX_N); return 1;}printf("Now input the (x_i,y_i),i=0,...,%d:\n",n);for(i=0;i<=n;i++)scanf("%lf%lf",&points[i].x,&points[i].y);printf("Now input the x value:");scanf("%lf",&x);for(i=0;i<=n;i++){diff[i]=0;tx=1;ty=1;for(j=0;j<=n;j++){if(i!=j){tx=tx*(x-points[j].x);ty=ty*(points[i].x-points[j].x);}}diff[i]=tx/ty;}for(i=0;i<=n;i++){tmp=points[i].y*diff[i];printf("%f",tmp);lagrange+=tmp;}printf("lagrange(%f)=%f\n",x,lagrange);return 0;}2、Newton 插值#include <stdio.h>#define MAX_N 20typedefstructtagPOINT{ double x;double y;} POINT;int main(){ int n;inti,j;POINT points[MAX_N+1];double diff[MAX_N+1];doublex,tmp,newton=0;printf("\nInput n value: ");scanf("%d",&n);if (n>MAX_N){printf("The input n is larger thenMAX_N,please redefine the MAX_N.\n");return 1;}if (n<=0){printf("Please input a number between 1 and %d.\n",MAX_N);return 1;}//输入被插值点(x_i,y_i)printf("Now input the (x_i,y_i),i=0,...,%d:\n",n);for (i=0;i<=n;i++)scanf("%lf%lf",&points[i].x,&points[i].y);printf("Now input the x value: ");scanf("%lf",&x);for (i=0;i<=n;i++)diff[i]=points[i].y;for (i=0;i<n;i++){for (j=n;j>i;j--){diff[j]=(diff[j]-diff[j-1])/(points[j].x-points[j-1-i].x);}}tmp=1;newton=diff[0];for(i=0;i<n;i++){tmp=tmp*(x-points[i].x);newton=newton+tmp*diff[i+1];}printf("newton(%f)=%f\n",x,newton);return 0;}四、运行结果1、Lagrange 插值1910年Larange插值计算得到的人口数:1965年Larange插值计算得到的人口数:2002年Larange插值计算得到的人口数:从插值计算得出的结果1910年的人口数是31872000人,1965年的人口数约为193081511人,2002年的人口数约为26138748,而1910年的实际人口数为91772000人,1960年的实际人口数为179323000人,1970年的人口数为203212000人,所以拉格朗日插值计算得出的结果只有1965年的人口数与实际值相差较近,而1910年和2002年的计算结果都与实际值相差较大,所以插值计算得到的数据准确性并不高。

牛顿插值法的C语言编程

牛顿插值法的C语言编程

Newton 插值Newton 插值函数 Newton 插值函数是用差商作为系数,对于01,,,n x x x …这1n +个点,其一般形式为:00100120101011()[][,]()[,,]()()[,,,]()()()n n n N x f x f x x x x f x x x x x x x f x x x x x x x x x −=+−+−−++−−−…………对于011,,,n x x x −…这n 个点,100100120101012()[][,]()[,,]()()[,,,]()()()n n n N x f x f x x x x f x x x x x x x f x x x x x x x x x −−=+−+−−++−−−…………差商的定义 若已知函数()f x 在点(0,1,2,,)i x i n =⋅⋅⋅处的函数值()i f x 。

则称:00[]()f x f x =为函数()f x 在点0x 的0阶差商;100110[][][,]f x f x f x x x x −=−为函数()f x 关于01,x x 的1阶差商;120101220[,][,][,,]f x x f x x f x x x x x −=−为函数()f x 过点012,,x x x 的2阶差商;依此类推,一般地称121012101210[,,,,][,,,,][,,,,,]k k k k k k k f x x x x f x x x x f x x x x x x x −−−−⋅⋅⋅−⋅⋅⋅⋅⋅⋅=−为函数()f x 关于01,,,k x x x ⋅⋅⋅的k 阶差商。

表1 差商表i x ()i f x 1阶差商 2阶差商3阶差商4阶差商0x 1x 2x 3x 4x……0()f x 1()f x 2()f x 3()f x 4()f x ……01[,]f x x12[,]f x x 23[,]f x x 34[,]f x x……012[,,]f x x x 123[,,]f x x x 234[,,]f x x x ……0123[,,,]f x x x x 1234[,,,]f x x x x ……01234[,,,,]f x x x x x……根据Newton 插值函数编写的C 语言编程根据Newton 插值函数并对照上面的差商表,可编写出Newton 插值法的C 语言程序如下: #include<stdio.h> #include<iostream.h> #include <stdlib.h>double NewtonInterpolation(double *x,double *y,int n,double xx,double *pyy) {double *f=(double *)malloc(n*sizeof(double));int i,k;for(i=1;i<=n-1;i++)f[i]=y[i];for(k=1;k<=n-1;k++)for(i=k;i<=n-1;i++){f[i]=(y[i]-y[i-1])/(x[i]-x[i-k]);if(i==n-1)for(i=k;i<n;i++)y[i]=f[i];}*pyy=y[n-1];for(i=n-2;i>=0;i--)*pyy=(*pyy)*(xx-x[i])+y[i];free(f);return 0;}void main(){int n=5;double x[5]={1.0,2.7,3.2,4.8,5.6},y[5]={14.2,17.8,22.0,38.3,51.7},xx=3,yy; NewtonInterpolation(x,y,n,xx,&yy);printf("%lf\n",yy);}。

Newton插值和三次样条插值的程序代码

Newton插值和三次样条插值的程序代码

(){}21()(11),5,10,20:12521()1,(0,1,2,,)()2,(0,1,2,,)()()235,20:1100(),i i i i n n k k k Newton f x x n x f x x i i n f x n x y i n Newton N x S x n x k y f x N =-≤≤=+=-+====-+=插值多项式和三次样条插值多项式。

已知对作、计算函数在点处的值;、求插值数据点的插值多项式和三次样条插值多项式;、对计算和相应的函数值()() (1,2,,99)4:()max()()max ()n k n k n k n k n k n k k k x S x k E N y N x E S y S x ==-=-和;、计算,;解释你所得到的结果。

Newton 插值的程序代码:%求Newton 插值多项式function new=newton(m,y)n=length(m);b=1;new=[zeros(1,n-1),y(1)];for i=2:n y(i:n)=(y(i:n)-y(i-1:n-1))./(m(i:n)-m(1:n-i+1)); %求差商b=conv(b,[1,-m(i-1)]);a=[zeros(1,n-length(b)),b];new=new+y(i)*a;end三次样条插值的程序代码:%求三次样条插值多项式function s=myspline(m,y)syms x;n=length(m);M=y;M(2:n)=(M(2:n)-M(1:n-1))./(m(2:n)-m(1:n-1));M(3:n)=(M(3:n)-M(2:n-1))./(m(3:n)-m(1:n-2));M=[0,M(3:n),0];h=m(2:n)-m(1:n-1);c=h(2:n-1)./(h(2:n-1)+h(1:n-2));a=1-c;b=2+zeros(1,n);M=6*M;c=[0,c];a=[a,0];M=chase(a,b,c,M); %用追赶法解方程组for i=1:n-1 %分区间计算插值多项式s(i)=M(i)*(m(i+1)-x)^3+M(i+1)*(x-m(i))^3+(6*y(i)-M(i)*h(i)^2)*(m(i+1)-x)+(6*y(i+1)-M(i+1)*h(i)^2 )*(x-m(i));s(i)=s(i)/(6*h(i));end主程序为:clear;clc;n=5; %指定n的值syms x %符号数xf=1./(1+25*x^2); %定义函数f(x)i=0:n;m=-1+2*i/n;y=subs(f,m); %计算函数值f(x)disp(['When n=',num2str(n),',the values of f(x) at the points of x(i) are:']);disp(y);%求Newton插值多项式new=newton(m,y);new=round(10000*new)/10000;nout=poly2sym(new);nout=vpa(nout,4); %整理多项式并输出disp(['When n=',num2str(n),',N(x)=']);disp(nout);%求三次样条插值多项式s=myspline(m,y);for i=1:n %分区间计算插值多项式a=sym2poly(s(i));b=poly2sym(a);b=vpa(b,4); %整理多项式并输出disp(['In the range of [',num2str(m(i)),',',num2str(m(i+1)),'],S(x)=']);disp(b);endk=0:100;xk=-1+0.02*k; %计算xk的值ll1=subs(f,xk);ll2=polyval(new,xk); %计算f(xk),N(xk)for i=1:length(xk)r=find(m(2:n+1)>=xk(i),1,'first');ll3(i)=subs(s(r),xk(i)); %计算S(xk)enden=max(abs(ll1-ll2));es=max(abs(ll1-ll3)); %计算E(Nn),E(Sn),并输出disp(['E(Nn)=',num2str(en)]);disp(['E(Sn)=',num2str(es)]);plot(xk,ll1,xk,ll2,'r',xk,ll3,'g');axis([-1,1,-0.1,1.1]);title('Interpolation');legend('Original function','Newton Interpolation','Spline')grid; %将原函数,Newton插值函数,三次样条插值函数画成图形。

c语言牛顿迭代法

c语言牛顿迭代法

c语言牛顿迭代法牛顿迭代法(Newton-Raphson法)是一种求解方程近似解的方法,它是利用泰勒级数展开函数在某点的值,然后用一阶泰勒展开式的根近似表示函数的零点,因此也被称为牛顿拉弗森法。

它可以高效地解决复杂的非线性方程组,是科学计算领域中最为常用和基础的方法之一。

牛顿迭代法的基本思想是:在第k次迭代时,求出曲线f(x)在点xk的一次导数斜率,以此确定x轴上的一个点xk+1,和该点处曲线的一次切线。

这条切线和x轴交点的横坐标就是极值点的估计值。

这个过程可以迭代多次,直到达到满足一定的误差精度或者迭代次数的要求。

C语言实现牛顿迭代法需要先定义一个函数,这个函数就是需要求解方程的函数。

定义完函数之后,需要实现牛顿迭代公式来求出下一次迭代的估计值,然后不断迭代。

具体实现过程如下:1. 定义函数f(x),即需要求解方程的函数。

2. 定义函数f_prime(x),即f(x)的一次导数。

3. 定义变量x和x_next,初始化它们的值。

4. 在循环中,首先计算f(x)和f_prime(x),然后计算下一个迭代点的估计值x_next = x - f(x) / f_prime(x)。

5. 如果x_next和x的差异满足预设的精度要求,则退出循环。

6. 否则,将x_next的值赋值给x,并重复执行第4步。

C语言实现牛顿迭代法的代码如下:#include <stdio.h>#include <math.h>定义函数f(x)double f(double x) {return x * x - 2;}定义函数f_prime(x)double f_prime(double x) {return 2 * x;}int main() {定义变量double x, x_next, epsilon;int iter;初始化变量x = 1.0;epsilon = 1e-6;iter = 0;迭代求解do {x_next = x - f(x) / f_prime(x);iter++;printf("Iteration %d: x = %lf\n", iter, x_next);x = x_next;} while (fabs(x_next - x) >= epsilon);输出结果printf("Final result: x = %lf\n", x);return 0;}在这个代码中,我们使用了do-while循环来不断执行迭代过程,直到达到预设的精度要求。

拉格朗日和牛顿插值法的C 方法实现(数值分析上机实验)

拉格朗日和牛顿插值法的C  方法实现(数值分析上机实验)

数值分析上机实验实验一一.上机题目:已知: 4 =2,9 =3,16 =4分别用二次Lagrange和Newton插值法求7 的近似值。

二.解题方法:1.lagrange方法:设x0=4,y0=2,x1=9,y1=3,x2=16,y2=4代入方程:(x1-X)(x2-X)/(x1-x0)(x2-x0)*y0+(x0-X)(x2-X)/(x0-x1)(x2-x1)*y1+(x1-X)(x0-X)/(x1-x2)(x0-x2)*y2令X=7代入方程得 Y=2.628572.Newton方法:设x0=4,y0=2,x1=9,y1=3,x2=16,y2=4建表4 29 3 0.216 4 0.14286 -0.00476f(x)=f(x0)+f[x0,x1](X-x0)+f[x0,x1,x2](X-x0)(X-x1)(X-x2)令X=7代入方程得Y=2.62857三.算法公式步骤:grange方法:通过公式写出算法并得出最后的值Y:for(b=0;b<m;b++)//完成公式f(Xn)外层嵌套循环f[b]=i//{double l=1;//保证每次跳出内层循环将L置1 不会将第一项的值带入下一项//for(a=0;a<m;a++)//完成公式f(Xn)内层嵌套循环f[a]=j//{if(a!=b)//完成定义i=1,i!=j//l=(f[a]-F)/(f[a]-f[b])*l;//完成(j-m)/(j-i)//la=l*g[b];//完成公式的F(X0)=f(X0)*Y0并累乘输出结果// }Y=la+Y;//累加x0y0+x1y1+...得最后结果//}2.Newton方法:先建表,通过二维数组的思想建表for(l=2;l<m+2;l++)//外层循环控制y阶数//{for(k=1;k<m+1;k++)//内层循环控制x个数//{a[k][l]=(a[k][l-1]-a[k-1][l-1])/(a[k][0]-a[k-l+1][0]);//完成f(x0,x1,...,xn)并存表//}}填表。

c语言牛顿法求方程的根

c语言牛顿法求方程的根

c语言牛顿法求方程的根牛顿法是一种迭代求根的方法,可以求解非线性方程的根。

下面是用C语言实现牛顿法求方程的根的示例代码:```c#include <stdio.h>#include <math.h>// 要求解的方程double f(double x) {return x*x - 2; // 求解x^2 - 2 = 0的根}// 方程的导数double df(double x) {return 2*x; // 方程的导数为2x}// 牛顿法求根函数double newton(double x0, double epsilon) {double x = x0; // 初始值double delta; // 误差do {double fx = f(x);double dfx = df(x);delta = fx / dfx; // 迭代公式x = x - delta; // 迭代计算新的x} while (fabs(delta) > epsilon); // 判断误差是否小于给定的精度return x;}int main() {double x0 = 1; // 初始值double epsilon = 0.0001; // 精度double root = newton(x0, epsilon);printf("方程的根为:%lf\n", root);return 0;}```以上代码中,`f`函数表示要求解的方程,`df`函数表示方程的导数,`newton`函数是用牛顿法求方程的根的函数,`main`函数是主函数,用于测试求解结果。

在代码中,初始值`x0`和精度`epsilon`可以根据实际情况进行调整。

运行代码得到的输出结果为方程的根值。

Newton法、一般迭代法Steffensen法、弦截法C语言代码

Newton法、一般迭代法Steffensen法、弦截法C语言代码

一、Newton法:#include<math.h>#include<stdio.h>double f(double x){return (3*x*x-exp(x));}double f1(double x){return (6*x-exp(x));}void main(){double x1=1,x;do{x=x1;x1=x-f(x)/f1(x);printf("x=%.9lf\n",x1);}while(fabs(x1-x)>0.000005);}说明:f 为原函数,f1为f的导函数,x1为初始值,通过"x=%.9lf“控制输入输出格式二、一般迭代法#include <stdio.h>#include <math.h>int main(){double x=1,x1;while(1){x1=pow(3*x+1,0.2);printf("x=%.6lf\n",x1);if(fabs(x1-x)<0.000005 )break;x=x1;}return 0;}说明:x1为初始值,x1=pow(3*x+1,0.2);为迭代格式,0.000005为允许误差,通过"x=%.6lf“控制输入输出格式三、Steffensen法:#include"stdio.h"#include"math.h"#define phi(x) pow(3*(x)+1,0.2);void main(){double x,x0,del,y,z;printf("x0="); scanf("%lf",&x0);printf("\ndel=:"); scanf("%lf",&del);while(1){y=phi(x0); z=phi(y);x=x0-(y-x0)*(y-x0)/(z-2*y+x0);printf("\n%.6lf",x);if(fabs(x-x0)<del) break;x0=x;}}说明:x0为初始值,pow(3*(x)+1,0.2);为φ(x)的格式,del为允许误差,通过"x=%.6lf“控制输入输出格式四、弦截法:#include<math.h>#include<stdio.h>double f(double x){ //计算f(x)的值return pow(2,x)+pow(3,x)-pow(4,x);}double point(double x1,double x2){//计算与x轴交点的x值printf("x=%.5f\n",(x1*f(x2)-x2*f(x1))/(f(x2)-f(x1)));return (x1*f(x2)-x2*f(x1))/(f(x2)-f(x1));}int main(){//输入两个数x1,x2double x1,x2,x;do{printf("输入两个数x1,x2:");scanf("%lf%lf",&x1,&x2);}while (f(x1)*f(x2)>= 0); // 当输入两个数大于0为真时,继续重新输入//关键循环步骤:do{x=point(x1,x2);//得到交点的值if(f(x)*f(x1)>0)x1=x;//新的x1elsex2=x;}while (fabs(f(x)) > 0.000005); }。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

Newton 插值
Newton 插值函数 Newton 插值函数是用差商作为系数,对于01,,,n x x x …这1n +个点,其一般形式为:
00100120101011()[][,]()[,,]()()[,,,]()()()
n n n N x f x f x x x x f x x x x x x x f x x x x x x x x x −=+−+−−++−−−…………对于011,,,n x x x −…这n 个点,
100100120101012()[][,]()[,,]()()[,,,]()()()
n n n N x f x f x x x x f x x x x x x x f x x x x x x x x x −−=+−+−−++−−−…………差商的定义 若已知函数()f x 在点(0,1,2,,)i x i n =⋅⋅⋅处的函数值()i f x 。

则称:
00[]()f x f x =为函数()f x 在点0x 的0阶差商;
100110
[][]
[,]f x f x f x x x x −=
−为函数()f x 关于01,x x 的1阶差商;
120101220
[,][,]
[,,]f x x f x x f x x x x x −=
−为函数()f x 过点012,,x x x 的2阶差商;
依此类推,一般地称
121012101210
[,,,,][,,,,]
[,,,,,]k k k k k k k f x x x x f x x x x f x x x x x x x −−−−⋅⋅⋅−⋅⋅⋅⋅⋅⋅=
−为函数()f x 关于01,,,k x x x ⋅⋅⋅的
k 阶差商。

表1 差商表
i x ()i f x 1阶差商 2阶差商
3阶差商
4阶差商
0x 1x 2x 3x 4x
……
0()f x 1()f x 2()f x 3()f x 4()
f x ……
01[,]f x x
12[,]f x x 23[,]f x x 34[,]f x x
……
012[,,]f x x x 123[,,]f x x x 234[,,]
f x x x ……
0123[,,,]f x x x x 1234[,,,]
f x x x x ……
01234[,,,,]f x x x x x
……
根据Newton 插值函数编写的C 语言编程
根据Newton 插值函数并对照上面的差商表,可编写出Newton 插值法的C 语言程序如下: #include<stdio.h> #include<iostream.h> #include <stdlib.h>
double NewtonInterpolation(double *x,double *y,int n,double xx,double *pyy) {
double *f=(double *)malloc(n*sizeof(double));
int i,k;
for(i=1;i<=n-1;i++)
f[i]=y[i];
for(k=1;k<=n-1;k++)
for(i=k;i<=n-1;i++)
{f[i]=(y[i]-y[i-1])/(x[i]-x[i-k]);
if(i==n-1)
for(i=k;i<n;i++)
y[i]=f[i];
}
*pyy=y[n-1];
for(i=n-2;i>=0;i--)
*pyy=(*pyy)*(xx-x[i])+y[i];
free(f);
return 0;
}
void main()
{
int n=5;
double x[5]={1.0,2.7,3.2,4.8,5.6},y[5]={14.2,17.8,22.0,38.3,51.7},xx=3,yy; NewtonInterpolation(x,y,n,xx,&yy);
printf("%lf\n",yy);
}。

相关文档
最新文档