数值分析拉格朗日实验题MATLAB程序,数值分析实验报告（拉格朗日插值牛顿插值最小二乘法）...

实验1

拉格朗日插值法

一、方法原理

n次拉格朗日插值多项式为：Ln(x)=y0l0(x)+y1l1(x)+y2l2(x)+…+ynln(x)

n=1时，称为线性插值，L1(x)=y0(x-x1)/(x0-x1)+

y1(x-x0)/(x1-x0)=y0+(y1-x0)(x-x0)/(x1-x0)

n=2时，称为二次插值或抛物线插值，精度相对高些

L2(x)=y0(x-x1)(x-x2)/(x0-x1)/(x0-x2)+y1(x-x0)(x-x2)/(x1-x0)/(x1-x2)+y2(x-x0)(x-x1)/(x2-x0)/(x2-x1)

二、主要思路

使用线性方程组求系数构造插值公式相对复杂，可改用构造方法来插值。

对节点xi(i=0,1,…,n)中任一点xk(0<=k<=n)作一n

次多项式lk(xk)，使它在该点上取值为1，而在其余点xi(i=0,1,…,k-1,k+1,…,n)上为0，则插值多项式为Ln(x)=y0l0(x)+y1l1(x)+y2l2(x)+…+ynln(x)

上式表明：n

个点xi(i=0,1,…,k-1,k+1,…,n)都是lk(x)的零点。可求得lk

三．计算方法及过程：1.输入节点的个数n

2.输入各个节点的横纵坐标

3.输入插值点

4.调用函数，返回z

函数语句与形参说明

程序源代码如下：

形参与函数类型

参数意义

int n

节点的个数

double x[n](double *x)

存放n个节点的值

double y[n](double *y)

存放n个节点相对应的函数值

double p

指定插值点的值

double fun()

函数返回一个双精度实型函数值，即插值点p处的近似函数值

#include

#include

using namespace std;

#define N

100

double fun(double *x,double *y, int n,double

p);

void main()

{int i,n;

cout<

cin>>n;

double x[N], y[N],p;

cout<

xiangliang x= "<

for(i=0;i>x[i];

cout<

xiangliang y= "<

for(i=0;i>y[i];

cout<

LagelangrichazhiJieDian p=

"<

cin>>p;

cout<

"<

system("pause") ;}

double fun(double x[],double y[], int n,double

p)

{double z=0,s=1.0;

int k=0,i=0;

double L[N];

while(k

{ if(k==0)

{

for(i=1;i

L[0]=s*y[0];

k=k+1;}

else

{s=1.0;

for(i=0;i<=k-1;i++)s=s*((p-x[i])/(x[k]-x[i]));

for(i=k+1;i

s=s*((p-x[i])/(x[k]-x[i]));

L[k]=s*y[k];

k++;}

}

for(i=0;i

return z;

}

四．运行结果测试：

五．实验分析

n=2时，为一次插值，即线性插值

n=3时，为二次插值，即抛物线插值

n=1，此时只有一个节点，插值点的值就是该节点的函数值

n<1时，结果都是返回0的；这里做了n=0和n=-7两种情况

3

常用的是线性插值和抛物线插值，显然，抛物线精度相对高些

n次插值多项式Ln(x)通常是次数为n的多项式，特殊情况可能次数小于n.例如：通过三点的二次插值多项式L2(x),如果三点共线，则y=L2(x)就是一条直线，而不是抛物线，这时L2(x)是一次式。

拟合曲线光顺性差

实验2

牛顿插值法

一、方法原理及基本思路

在拉格朗日插值方法中，若增加一个节点数据，其插值的多项式需重新计算。现构造一个插值多项式Nn(x)，只需对Nn-1(x)作简单修正(如增加某项)即可得到，这样计算方便。

利用牛顿插值公式，当增加一个节点时，只需在后面多计算一项，而前面的计算仍有用；另一方面Nn(x)的各项系数恰好又是各阶差商，而各阶差商可用差商公式来计算。

由线性代数知，对任何一个不高n次的多项式P(x)=b0＋b1x＋b2x2＋…＋bnxn (幂基) ①

也可将其写成P(x)=a0＋a1(x－x0)＋a2(x－x0)

(x－x1)＋…＋an(x－x0)

…(x－xn-1)

其中ai为系数，xi为给定节点，可由①求出ai 一般情况下，牛顿插值多项式Nn(x)可写成：

Nn(x)=

a0＋a1(x－x0)＋a2(x－x0)

(x－x1)＋…＋an(x－x0)

…(x－xn-1))

只需求出系数ai，即可得到插值多项式。

二、计算方法及过程

1.先后输入节点个数n和节点的横纵坐标，插值点的横坐标，最后输入精度e

2. 用do-while循环语句得到跳出循环时k的值

3.将k值与n-1进行比较，若在达到精度时k

函数语句与形参说明

形参与函数类型

参数意义

int n

节点的个数

float x[MAX]

存放n个节点的值(从小到大)

Float y[MAX];

存放n个节点相对应的函数值

float x0,y0

指定插值点的横纵坐标

float e

精度

程序源代码如下：

#include

#include

using namespace std;

#define MAX

100

void main()

{ float x[MAX],y[MAX];

float

x0,y0,e,N1;

float

N0=0;

int i,k,n;

cout<

cin>>n;

cout<

for(i=0;i

cin>>x[i];

cout<

for(i=0;i

cin>>y[i];

cout<

cin>>x0;

cout<

cin>>e;

y0=1; N1=y[0]; k=0;

do { k++;

N0=N1;

y0=y0*(x0-x[k-1]);

for(i=0;i

y[k]=(y[k]-y[i])/(x[k]-x[i]);

N1=N0+y0*y[k];

} while (fabs(N1-N0) > e

&& k

if (k==(n-1))

cout<

if

(k

cout<

system("pause");

}

三．运行结果测试：

题一：已知f(x)=sh(x)的函数表如下：计算f(0.23)的近似值

xi

0

0.20

0.30

0.50

Sh(xi)

0

0.20134

0.30452

0.52110

题二：已知根号100等于10，根号121等于11，根号144等于12； 用二次插值和三次插值计算根号115的近似值。

四．实验分析

显然，题二中当n=2,n=3时，k的值没有在n-1前满足各自的精度要求，所以计算失败

从题一中可以看出，插值点的个数、精度、插值点的选择都会影响实验的结果；

我们通常会选择与插值点最接近的节点，可以提高精度；

在可以计算出结果的情况下，插值点越多，结果越精确。

实验3

最小二乘法求拟合曲线

一、方法原理及思路

已知数据对(xj, yj)(j=1,2,…,n)，求多项式

为最小，这就是一个最小二乘问题。

二、计算方法及过程：

1. 先后输入数据的对数m,近似多项式的最高次数n及各个数据的横纵坐标

2. 用for循环语句得到线性方程组Ga=d

3. 用高斯列主元消去法得到的解即为近似多项式的系数ai(i=0,1,…,n)

函数语句与形参说明

形参与函数类型

参数意义

int m

数据的对数,即坐标的个数

int n

近似多项式的最高次数

float x[MAX]

存放m个横坐标的值(从小到大)

Float y[MAX];

存放m个纵坐标的值(与相应的横坐标对应)

Float a[ MAX +1],

近似多项式的系数

float G[ MAX +1][ MAX +1]

法方程的系数矩阵，且为对称的，正定的

float s[2 * MAX +1]

系数矩阵G中的元素

Float b[ MAX +1];

法方程中等号右边的列向量

程序原代码如下

#include

#include

#define MAX 100

using namespace std;

int

main()

{ float

t;

int

i,j,k,r,m,n;

cout<

cin>>m;

cout<

cin>>n;

float y[

MAX],x[ MAX ],x1[ MAX ],a[ MAX +1],G[ MAX +1][ MAX +1], s[2* MAX

+1],b[ MAX +1];

cout<

for

(i=0;i>x[i];

cout<

for

(i=0;i>y[i];

s[0]=m;

b[0]=0.0;

for(j=0;j

{b[0]=b[0]+y[j];

x1[j]=1.0;

}

for

(i=1;i<=n;i++)

{b[i]=0.0; s[2*i-1]=0.0;s[2*i]=0.0;

for(j=0;j

{x1[j]=x1[j]*x[j];

b[i]=b[i]+x1[j]*y[j];

s[2*i-1]=s[2*i-1]+

x1[j]*x1[j]/x[j];

s[2*i]=s[2*i]+x1[j]*x1[j];}

}

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

{ for(j=0;j<=n;j++)

G[i][j]=s[i+j]; }

for(k=0;k

{ t=abs(G[k][k]);

r=k;

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

{

if(t

{ t=G[i][k];

r=i; }

}

if(r>k)//行交换

{ for(j=k;j<=n;j++)

{ t=G[k][j];

G[k][j]=G[r][j];

G[r][j]=t; }

t=b[k]; b[k]=b[r];

b[r]=t; }

for(i=k+1;i<=n;i++)//消元

{

t=G[i][k]/G[k][k];

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

G[i][j]=G[i][j]-t*G[k][j];

b[i]=b[i]-t*b[k];

}

}

a[n]=b[n]/G[n][n];//回代

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

{ i=n-k; a[i]=b[i];

for(j=i+1;j<=n;j++)a[i]=a[i]-G[i][j]*a[j];

a[i]=a[i]/G[i][i]; }

cout<

//输出答案

for(i=0;i<=n;i++)cout<

system("pause");

}

三：运行结果测试：

四．实验分析

由于曲线拟合的最小二乘法一般是经过描点，确定其近似多项式的形式，但由于给定的点有误差，所以拟合曲线的数学模型并不是一开始就能选的好的，往往要通过分析确定若干模型后，再经过实际计算，才能选到较好的模型。

对我的博文感兴趣的话，关注我哦~~~~