Dwcsoic数值分析实验报告1.docx

1、Dwcsoic数值分析实验报告1七夕，古今诗人惯咏星月与悲情。吾生虽晚，世态炎凉却已看透矣。情也成空，且作“挥手袖底风”罢。是夜，窗外风雨如晦，吾独坐陋室，听一曲尘缘，合成诗韵一首，觉放诸古今，亦独有风韵也。乃书于纸上。毕而卧。凄然入梦。乙酉年七月初七。-啸之记。 For personal use only in study and research; not for commercial use数值分析实 验 报 告 册姓名: 学号: 专业: 年级: 计算机科学学院计算机应用教研室 2008 年 春季 学期实验七18 实验一一、课题名称非线性方程数值解法二、目的和意义学会常用的插值方法，求函

2、数的近似表达式，以解决其它实际问题；明确插值多项式和分段插值多项式各自的优缺点；熟悉插值方法的程序编制；如果绘出插值函数的曲线，观察其光滑性。三、计算公式Lagrange插值公式： 牛顿插值公式：四、结构程序设计 程序设计:#includemath.hfloat f(float x)return(x*x*x-1)/3); /*牛顿迭代函数*/main()float x1,x2,eps,d;int k=0;clrscr();printf(n input x1=); /*输入迭代初值*/scanf(%f,&x1);printf(n input eps=); /*输入求解精度eps*/scanf(%

3、f,&eps);do k+; x1=x2; x2=f(x1); printf(n %d %fn,k,x2); while(fabs(x2-x1)=eps); printf(the root of f(x)=0 is x=%f,k=%dn,x2,k); /*输出x和迭代次数k*/ getch(); 五、结果讨论和分析 计算结果分析：将六种迭代格式分别代入程序试验：（1）第一种格式：无论何值都无法求出，即发散（2）第二种格式：初值为任意的x（x20），精度为0.0001 X=1.879372，k=10 其他值为发散。（4）第四种格式：初值为任意值，精度为0.00001 X=-0.347296,k=

4、5（5）第五种格式：初值为任意值,精度为0.00001 X=-0.347296,k=4（6）第六种格式：初值为任意值，精度为0.00001 X=-0.347296,k=4 由此可知不同的初值对公式的计算有影响，当初值不满足函数的收敛条件时，无法计算结果，函数发散。 精度的大小不同也使迭代函数迭代的次数不同，从而影响xn的近似程度。实验二一、课题名称解线性方程组的直接方法二、目的和意义 掌握线性方程组直接接法的基本思想；了解不同数值方法解线性方程组的原理、实现条件、使用范围、计算公式；培养编程与上机调试能力。三、计算公式 消去法 设a(k)kk=0,对k=1，2，n-1计算 mik=a(k)ik

5、/a(k)kk a(k+1)ij=a(k)ij-mika(k)kj i,j=k+1,k+2,n b(k+1)i=b(k)i-mikb(k)kn xn=b(n)n/a(n)nnj=i+1 xi=(b(i)i-a(i)ijxj)/a(i)ii i=n-1,n-2,1 平方根法 追赶法 lij=(aii-l2ik)1/2 Ly=f lji=(aji-ljklik)/lii j=i+1,i+2,n Ux=y y1=f1/l1 y2=(fi-aiyi-1)/li i=2,3,n 四、结构程序设计用追赶法求解线性方程组#includestdio.hmain() FILE*f; double a15,b15

6、,c15,d15; double t; int i,n; f=fopen(zgf.dat,r); fscanf(f,%d,&n); fscanf(f,%lf%lf%lf,&b1,&c1,&d1); for(i=2;i=n-1;i+) fscanf(f,%lf%lf%lf%lf,&ai,&bi,&ci,&di); fscanf(f,%lf%lf%lf,&an,&bn,&dn); fclose(f); c1=c1/b1; d1=d1/b1; for(i=2;i=1;i-)di=di-ci*di+1; printf(n*n); for(i=1;i=n;i+) printf(d%2d=%lfn,i,d

7、i);五、结果讨论和分析 此方法通过有限步算术运算求出 精确解，但实际计算由于舍入误差的影响，只能求出近似解。实验三一、课题名称解线性方程组的迭代法二、目的和意义 了解各迭代法的基本原理和特点，判断雅克比迭代、高斯-塞德尔迭代对任意初始向量的收敛性，完成雅克比迭代、高斯-塞德尔迭代算法的程序实现三、计算公式 雅可比 xi(k+1)=1/aii(bi-aijxj(k) 高斯-塞德尔xi(k+1)=1/aii(bi-aijxj(k+1)-aijxj(k) 超松弛迭代 xi(k+1)=（1-w）xi(k)+w*(bi-aijxj(k+1)-aijxj(k) /aii四、结构程序设计高斯-塞德尔法：#

8、includemath.h#define M 8#define N 9main() double aMN=4,2,-4,0,2,4,0,0,0, 2,2,-1,-2,1,3,2,0,-6, -4,-1,14,1,-8,-3,5,6,20, 0,-2,1,6,-1,-4,-3,3,23, 2,1,-8,-1,22,4,-10,-3,9, 4,3,-3,-4,4,11,1,-4,-22, 0,2,5,-3,-10,1,14,2,-15, 0,0,6,3,-3,-4,2,19,45; double xM=0,0,0,0,0,0,0,0; double r,t,q,eps=0.0001; int k,

9、i,j,T=100; for(i=0;iM;i+) for(j=0;jN;j+) printf(%1f,aij); printf(n); for(k=0;kT;k+) r=0; for(i=0;iM;i+) t=xi; q=0; for(j=0;jr)r=fabs(xi-t); if(reps)break; printf(nk=%d,k); for(i=0;iM;i+) printf(nx%d=%lf,i,xi); if(k=T)printf(nNo); else for(i=0;iM;i+) printf(x(%d)=%15.7fn,i+1,xi); 五、结果讨论和分析与直接法相比，迭代法适

10、用于稀疏矩阵的线性方程组实验四一、课题名称函数插值方法二、目的和意义 了解多项式差值公式的存在唯一性条件及其余项表达式的推导，了解拉格朗日插值多项式的构造、计算及其基函数的特点，牛顿插值多项式的构造与应用，差商、差分的计算及基本性质。三、计算公式=，i=0,1,2nP(X)=p(x)+ (p(x)的初值是0)。四、结构程序设计#includestdio.h#includemath.h#includestring.h#includeconio.h#includestdlib.h#define n 5double x1=0.4,0.55,0.65,0.80,0.95,1.05;double y1=

11、0.41075,0.57815,0.69675,0.90,1.00,1.25382;main() double Lag(double x1,double y1,float t); int m,k;float x,y;float X;double z; printf(nthe number of the interpolation points is m:); scanf(%d,&m); for(k=1;k=m;k+) printf(ninputX%d=,k); scanf(%f,&X); z=Lag(x1,y1,X); printf(P(%f)=%fn,X,z); getch(); retur

12、n(0);double Lag(double x,double y,float X) int i,j; double L,P; P=0.0; for(i=0;i=n;i+) L=1.0; for(j=0;j=n;j+) if(j!=i) L=L*(X-xj)/(xi-xj); P=P+yi*L; return (P);五、结果讨论和分析实验五一、课题名称曲线拟合的最小二乘法二、的和意义掌握曲线拟合的最小二乘法；了解最小二乘法亦可以用于解超定线性方程组；探索拟合函数的选择与拟合精度间的关系。三、计算公式e 22=2i=(xi)-f(xi)2四、结构程序设计#include stdio.h#inc

13、lude math.h#define num 10float neiji(float bnum,float cnum) int p;float nj=0;for(p=1;pnum;p+)nj+=cp*bp;return nj;float snum,xnum,fainumnum,afanum;float beidanum,anum,xfainum,ydnum,max,pcpfh;void main() int i,j,k,n,index,flag;char conti;conti= ;printf(请输入已知点的个数n=n);scanf(%d,&n);printf(请输入x和y:);for(i=

14、1;i=n;i+) printf(x%d=,i);scanf(%f,&xi);printf(y%d=,i);scanf(%f,&yi);while(conti= ) printf(请输入拟和次数=);scanf(%d,&index);pcpfh=0;afa1=0;a0=0;for(i=1;i=n;i+) afa1+=xi;a0+=ydi;fai0i=1;afa1=afa1/n;a0=a0/n;for(i=1;i=n;i+)fai1i=xi-afa1;a1=neiji(fai1,yd)/neiji(fai1,fai1);for(k=1;kindex;k+) for(i=1;i=n;i+)xfai

15、i=xi*faiki;afak+1=neiji(faik,xfai)/neiji(faik,faik);beidak=neiji(faik,faik)/neiji(faik-1,faik-1);for(j=1;j=n;j+)faik+1j=(xj-afak+1)*faikj-beidak*faik-1j;ak+1=neiji(faik+1,yd)/neiji(faik+1,faik+1);printf(%d次拟和结果为n,index);for(i=0;i=index;i+)printf(a%d=%fn,i,ai);for(i=1;i=index;i+)printf(afa%d=%fn,i,af

16、ai);for(i=1;iindex;i+)printf(beida%d=%fn,i,beidai);for(i=1;i=n;i+) for(k=0;k=index;k+)si+=ak*faiki;ydi=fabs(ydi-si);pcpfh+=ydi*ydi;si=0;max=0;for(i=1;imax)max=ydi;flag=i;printf(当x=%f时,偏差最大=%f,偏差平方和为%fn,xflag,max,pcpfh);printf(继续拟和请按space,按其他键退出);conti=getchar();conti=getchar();五、结果讨论和分析请输入已知点的个数n=10

17、请输入x和y：x1=0y1=0请输入拟合次数=55次拟合结果为当x=0.000000时，偏差最大=6706185.000000，偏差平方和为6126336.000000.实验六一、课题名称数值积分与数值微分二、目的和意义深刻认识数值积分法的意义；明确数值积分精度与步长的关系；根据定积分的计算方法，可以考虑二重积分的计算问题。三、计算公式=1/6*f(a)+4f(xk+1/2)+2f(xk)+f(b) Rn=64/63*c2n-1/63*cn四、结构程序设计Romberg算法:#includestdio.h#includemath.h#includeconio.hfloat f(float x)

18、return(exp(x)/(4+x*x);main() float a=0,b=1; float h=b-a,T1,T2,s,x; int i;T1=h/2*(f(a)+f(b); for(i=0;i3;i+) printf(h=%f,T=%fn,h,T1); s=0; x=a+h/2; while(xb) s+=f(x); x+=h; T2=T1/2+h/2*s; h/=2;T1=T2; printf(h=%f,T=%fn,h,T1); return;Simpson算法：#includestdio.h#includemath.h#includeconio.h#define Max_M 20

19、float f(float x)return(sin(x)/x);float Simpson(float a,float b,int n) int k;float x,s1,s2,h=(b-a)/n; x=a+h/2; s1=f(x);s2=0; for(k=1;kn;k+) s1=s1+f(a+k*h+h/a); s2=s2+f(a+k*h); s2=h*(f(a)+4*s1+2*s2+f(b)/6; return(s2);main() int i,n;float a1,b1,s=0; printf(n Input the begin:); scanf(%f,&a1); printf(n I

20、nput the end:); scanf(%f,&b1); do printf(n input n divde valuedivide(%f,%f):,a1,b1); scanf(%d,&b1); while(nMax_M); s=Simpson(a1,b1,n); printf(solve is:%f,s); getch(); return(s);_五、结果讨论和分析用Romberg法得出结果为：用Simpson法得出结果为：可见复化公式要先估计出步长，步长的大小将影响计算结果和精度 实验七一、课题名称常微分方程的数值解法二、目的和意义熟悉各种初值的问题的算法，编出算法程序；明确各种算法的

21、精度寓所选步长有密切关系；通过计算更加了解各种算法的优越性；三、计算公式k1=f(xi,yi)k2=f(xi+1/2+1/2,yi+h/2*k1)k3=f(xi+1/2,yi+h/2*k2)k4=f(xi+1,yi+h*k3)Yi+1=yi+h*(k1,2*k2+2*k3+k4)/64、结构程序设计 Rung-kutta法：#includemath.h#includestring.h#includestdio.h#includeconio.hfloat f(float x,float y) float y1; y1=y-2*x/y; return y1;float Runge_Kutta(fl

22、oat x,float y,float h)float k1,k2,k3,k4;k1=f(x,y);k2=f(x+h/2,y+h*k1/2);k3=f(x+h/2,y+h*k2/2);k4=f(x+h,y+h*k3);return(y+h*(k1+2*k2+2*k3+k4)/6);main() int i=0; float x,y,h,b; clrscr(); printf(n Input begin x0:); scanf(%f,&x); printf(n Input begin y0:); scanf(%f,&y); printf(n Input step h:); scanf(%f,&h

23、); printf(n Input end b:); scanf(%f,&b); printf(n x0=%10f y0=%10fn,x,y); do y=Runge_Kutta(x,y,h); x=x+h;i+; printf(x%d=%10f y%d=%10fn,i,x,i,y); while(xb); getch(); return(y); 5、结果讨论和分析 可见步长的选取会影响节点处数值解的误差 仅供个人用于学习、研究；不得用于商业用途。For personal use only in study and research; not for commercial use.Nur fr den persnlichen fr Studien, Forschung, zu kommerziellen Zwecken verwendet werden.Pour l tude et la recherche uniquement des fins personnelles; pas des fins commerciales. , , . 以下无正文