计算机上机报告.docx

计算机上机报告.docx

《计算机上机报告.docx》由会员分享，可在线阅读，更多相关《计算机上机报告.docx（23页珍藏版）》请在冰豆网上搜索。

计算机上机报告

计算方法上机实验报告

上课时间：

2014-2015学年秋学期，6～14周

一．拉格朗日插值------------------------------------------------------1

二．牛顿插值------------------------------------------------------------3

三．改进欧拉法---------------------------------------------------------5

四．四阶龙格-库塔-----------------------------------------------------7

五．牛顿迭代------------------------------------------------------------9

6．复化Simpson公式------------------------------------------------11

七．Romberg算法------------------------------------------------------14

八.Seidel迭代法------------------------------------------------------17

九.Gauss列主元消去法----------------------------------------------20

一.拉格朗日插值

1.程序代码

#include

voidLagrange（）

{

inti=0;

doublea[10],b[10],L,L1,L2,L3,L4,x;

cout<<"x=";

for（i=0;i<4;i++）

{

cin>>a[i];

}

cout<<"y=";

for（i=0;i<4;i++）

{

cin>>b[i];

}

cout<<"x=";

cin>>x;

L1=（x-a[1]）*（x-a[2]）*（x-a[3]）*b[0]/（a[0]-a[1]）/（a[0]-a[2]）/（a[0]-a[3]）;

L2=（x-a[0]）*（x-a[2]）*（x-a[3]）*b[1]/（a[1]-a[0]）/（a[1]-a[2]）/（a[1]-a[3]）;

L3=（x-a[0]）*（x-a[1]）*（x-a[3]）*b[2]/（a[2]-a[0]）/（a[2]-a[1]）/（a[2]-a[3]）;

L4=（x-a[0]）*（x-a[1]）*（x-a[2]）*b[3]/（a[3]-a[0]）/（a[3]-a[1]）/（a[3]-a[2]）;

L=L1+L2+L3+L4;

cout<<"L="<

}

voidmain（）

{

Lagrange（）;

cout<

}

2.例子

3.运行结果

二．牛顿插值

1.程序代码

#include

#include

voidmain（）

{

intn,i,j;

doubleA[50][50],*x,*y;

cout<<"请输入插值节点数：

";

cin>>n;

x=newdouble[n];

y=newdouble[n];

cout<<"请输入这"<

"<

for（i=0;i<=n-1;i++）

cin>>x[i]>>y[i];

doubleK=1,xx,N=0,P;

for（i=0;i<=n-1;i++）

{

A[i][0]=x[i];

A[i][1]=y[i];

}

for（j=2;j<=n;j++）

for（i=1;i<=n-1;i++）

{

A[i][j]=（A[i][j-1]-A[i-1][j-1]）/（A[i][0]-A[i-j+1][0]）;

}

for（i=0;i<=n-1;i++）

cout<<"输出第"<

"<

cout<<"请输入预求值x=";

cin>>xx;

for（i=0;i

{

K*=xx-x[i];

N+=A[i+1][i+2]*K;

P=A[0][1]+N;

}

cout<<"Newton插值结果为：

y="<

getch（）;

}

2.例子

3.运行结果

三．改进欧拉法

1.程序代码

#include

#include

doublefun（doublex,doubley）

{

return（-0.9*y/（1+2*x））;

}

voidmain（）

{

doublea,b,*y,h,*x,yp,yc;

intn,k;

cout<<"常微分方程为y'=-0.9*y/（1+2*x）"<

cout<<"其中0<=x<=1"<

cout<<"初值为y（0）=1"<

cout<<"请输入计算区间（a,b）：

";

cin>>a>>b;

cout<<"请输入步长h：

";

cin>>h;

cout<<"请输入计算次数：

";

cin>>n;

y=newdouble[n];

x=newdouble[n];

cout<<"请输入初值y（0）=";

cin>>y[0];

x[0]=a;

for（k=0;k<=n;k++）

{

yp=y[k]+h*fun（x[k],y[k]）;

yc=y[k]+h*fun（x[k]+h,yp）;

y[k+1]=0.5*（yp+yc）;

x[k+1]=x[k]+h;

}

cout<<"迭代结果为：

"<

for（k=0;k<=n;k++）

cout<<"x"<

getch（）;

}

2.例子

3.运行结果

四．四阶龙格-库塔

1.程序代码

#include

#include

doublefun（doublex,doubley）

{

return（x-y）;

}

voidmain（）

{

doublea,b,*y,h,x,k1,k2,k3,k4;

intn,k;

cout<<"常微分方程为y'=x-y"<

cout<<"其中0<=x<=1"<

cout<<"初值为y（0）=0"<

cout<<"请输入计算区间（a,b）：

";

cin>>a>>b;

cout<<"请输入步长h：

";

cin>>h;

cout<<"请输入计算次数：

";

cin>>n;

y=newdouble[n];

cout<<"请输入初值y（0）：

";

cin>>y[0];

x=a;

cout<<"Runge-Kutta迭代法结果为：

"<

cout<<"x0="<

for（k=0;k<=n-1;k++）

{

k1=fun（x,y[k]）;

k2=fun（x+h/2,y[k]+k1*h/2）;

k3=fun（x+h/2,y[k]+k2*h/2）;

k4=fun（x+h,y[k]+k3*h）;

y[k+1]=y[k]+（h/6）*（k1+2*（k2+k3）+k4）;

cout<<"x"<

x=x+h;

}

getch（）;

}

2.例子

3.运行结果

五．牛顿迭代法

1.程序代码（C++代码）

#include

#include

usingnamespacestd;

doublenewtondiedai（doublea,doubleb,doublec,doubled,doublex）;

intmain（）

{

doublea,b,c,d;

doublex=1.5;

cout<<"请依次输入方程四个系数：

";

cin>>a>>b>>c>>d;

x=newtondiedai（a,b,c,d,x）;

cout<

return0;

}

doublenewtondiedai（doublea,doubleb,doublec,doubled,doublex）

{

while（abs（a*x*x*x+b*x*x+c*x+d）>0.000001）

{

x=x-（a*x*x*x+b*x*x+c*x+d）/（3*a*x*x+2*b*x+c）;

}

returnx;

}

2.例子

3.运行结果

六．复化Simpson公式

1.程序代码（C++代码）

#include

#include

doublefunction1（doublex）//被积函数

{

doubles;

s=x/（4+x*x）;

returns;

}

doublefunction2（doublex）//被积函数

{

doubles;

s=sqrt（x）;

returns;

}

doubleReiterationOfSimpson（doublea,doubleb,doublen,doublef（doublex））//复化Simpson公式

{

doubleh,fa,fb,xk,xj;

h=（b-a）/n;

fa=f（a）;

fb=f（b）;

doubles1=0.0;

doubles2=0.0;

for（intk=1;k

{

xk=a+k*h;

s1=s1+f（xk）;

}

for（intj=0;j

{

xj=a+（j+0.5）*h;

s2=s2+f（xj）;

}

doublesn;//和

sn=h/6*（fa+fb+2*s1+4*s2）;//复化Simpson公式

returnsn;

}

main（）

{doublea,b,Result,n;

cout<<"请输入积分下限:

"<

cin>>a;

cout<<"请输入积分上限:

"<

cin>>b;

cout<<"请输入分割区间数n:

"<

cin>>n;

cout<<"复化Simpson公式计算结果:

";

Result=ReiterationOfSimpson（a,b,n,function1）;

cout<

}

2.例子（

）

3.运行结果

七．Romberg算法

1.程序代码（C++代码）

#include

#include

usingnamespacestd;

#definef（x）（4/（1+x*x））

#defineepsilon0.0001

#defineMAXREPT10

doubleRomberg（doubleaa,doublebb）

{intm,n;

doubleh,x;

doubles,q;

doubleep;

double*y=newdouble[MAXREPT];

doublep;

h=bb-aa;

y[0]=h*（f（aa）+f（bb））/2.0;

m=1;

n=1;

ep=epsilon+1.0;

while（（ep>=epsilon）&&（m

{p=0.0;

for（inti=0;i

{x=aa+（i+0.5）*h;

p=p+f（x）;

}

p=（y[0]+h*p）/2.0;

s=1.0;

for（intk=1;k<=m;k++）

{s=4.0*s;

q=（s*p-y[k-1]）/（s-1.0）;

y[k-1]=p;

p=q;

}

p=fabs（q-y[m-1]）;

m=m+1;

y[m-1]=q;

n=n+n;h=h/2.0;

}

return（q）;

}

intmain（）

{doublea,b;

cout<<"Romberg积分,请输入积分范围a,b:

"<

cin>>a>>b;

cout<<"积分结果:

"<

system（"pause"）;

return0;

}

2.例子（

）

3.运行结果

八．Seidel迭代法

1.程序代码（C++代码）

#include

#include

#definemax100

#defineEPS1e-6

floata[3][3]={{10,-1,-2},{-1,10,-2},{-1,-1,5}};

floatb[3]={7.2,8.3,4.2};

floatx[3]={0,0,0};

floaty[3];

floatS（intm）

{intn;

floatS=0;

floaty;

for（n=0;n<3;n++）

{if（m==n）{}

else

{S+=a[m][n]*x[n];

}

}y=（b[m]-S）/a[m][m];

returny;

}

voidmain（）

{inti;

intF,T=1,k=1;

do

{F=0;

if（T）

{T=0;}

else

{k++;

for（i=0;i<3;i++）

{x[i]=y[i];

}

}

for（i=0;i<3;i++）

{y[i]=S（i）;

}

for（i=0;i<3;i++）

{if（fabs（x[i]-y[i]）>EPS）

{F=1;

}

}printf（"%d\n",k）;

}while（（（F==1）&&k!

=max））;

printf（"迭代次数:

%d\n",k）;

for（i=0;i<3;i++）

{printf（"x[%d]=%1.6f\n",i,y[i]）;

}

}

2.例子（

）

3.运行结果

九．Gauss列主元消去法

1.程序代码（C++代码）

#include

#include

#include

#definemax_dimension20

intn;

staticfloata[max_dimension][max_dimension];

staticfloatb[max_dimension];

staticfloatx[max_dimension];

voidmain（）

{inti;

intj;

intd;

introw;

floattemp;

floatknown_items;

floatl[max_dimension][max_dimension];

printf（"请输入阶数:

"）;

scanf（"%d",&n）;

printf（"\n"）;

printf（"请输入系数矩阵的值:

"）;

printf（"\n"）;

for（i=0;i

{printf（"输入第%d行的值:

",i+1）;

for（j=0;j

{

scanf（"%f",&a[i][j]）;

}

printf（"\n"）;

}

printf（"请输入常数项的值:

"）;

for（i=0;i

scanf（"%f",&b[i]）;

printf（"\n"）;

for（d=0;d

{row=d;

for（i=d+1;i

{if（fabs（a[i][d]）>fabs（a[row][d]））

row=i;

}

if（row!

=d）

{for（j=d;j

{temp=a[row][j];

a[row][j]=a[d][j];

a[d][j]=temp;

}

temp=b[row];

b[row]=b[d];

b[d]=temp;

}

for（i=d+1;i

{l[i][d]=-a[i][d]/a[d][d];

for（j=d;j

{a[i][j]=a[i][j]+a[d][j]*l[i][d];

}

b[i]=b[i]+b[d]*l[i][d];

}

}

printf（"\n"）;

for（i=n-1;i>-1;i--）

{known_items=0;

for（j=1;j

{known_items=known_items+a[i][i+j]*x[i+j];

}

x[i]=（b[i]-known_items）/a[i][i];

}

printf（"方程组的根为:

\n\n"）;

for（i=0;i

printf（"%.5f",x[i]）;

printf（"\n\n"）;

getch（）;

}

2.例子

3.运行结果