数值计算方法编程作业C语言版

1、实用文案1:第二章(1)二分法求解非线性方程：#in clude<stdio.h>#in elude <math.h>#defi ne f(x) (x*x-1)*x-1)void mai n() float a,b,x,eps;int k=0;prin tf("i ntput epsn");/*容许误差 */sca nf("%f", &eps);prin tf("a,b=n ”);for(;)scanf("%f, %f',&a，&b);if(f(a)*f(b)>=0) /*

2、判断是否符合二分法使用的条件*/printf(”二分法不可使用，请重新输入:n");else break;dox=(a+b)/2;k+;if(f(a)*f(x)<0) /* 如果f(a)*f(x)<0，则根在区间的左半部分*/b=x;else if(f(a)*f(x)>0)/*否则根在区间的右半部分*/a=x;else break;while(fabs(b-a)>eps);/*判断是否达到精度要求，若没有达到，继续循环*/x=(a+b)/2; /*取最后的小区间中点作为根的近似值*/printf("n The root is x=%f, k=%dn&

3、quot;,x,k);运行结果：in tput eps0.00001a,b=2,-5The root is x=1.324721, k=20Press any key to con ti nue总结：本题关键在于两个端点的取值和误差的判断，此程序较容易。二分法收敛速度较快, 但缺点是只能求解单根。(2)牛顿法求解非线性方程：#i nclude <stdio.h>#in clude <math.h>float f(float x) /* 定义函数 f(x) */return(-3*x+4)*x-5)*x+6;float f1(float x) /* 定义函数 f(x)的导

4、数 */ return (-9*x+8)*x-5;void mai n() float eps,x0,x1=1.0;prin tf("i nput eps:n");scan f("%f", &eps); /*输入容许误差*/do x0=x1;/*准备下一次迭代的初值*/x1=x0-f(x0)/f1(x0);/* 牛顿迭代 */while(fabs(x1-x0)>eps); /*当满足精度，输出近似根*/prin tf("x=%fn",x1);程序运行结果：x=1.265328总结：关键是牛顿迭代的应用，程序中最大缺点是函

5、数及其导数已唯一给出确定不可求的随意函数的根，牛顿法比二分法收敛快，可以求重根。2:第三章(1)列主元素消去法求解线性方程：#in clude<iostream>#in clude<cmath>#defi ne N 20using n amespace std;void load();float aNN;int m;int mai n()int i,j;int c,k,n ,p,r;float xN,lNN,s,d;cout<<"下面请输入未知数的个数m="cin»m;cout<<e ndl;cout<<

6、"请按顺序输入增广矩阵a:"<<e ndl;load();for(i=0;i<m;i+)for( j=i;j<m;j+)c=(fabs(aji)>fabs(aii)?j:i;/* 找列最大元素 */for(n=0;n< m+1; n+)s=ain; ain=acn; acn=s;/* 将列最大数防在对角线上*/for(p=0;p<m+1;p+)cout<<aip<<"t"cout<<e ndl;for(k=i+1;k<m;k+) lki=aki/aii;for(r=i;r

7、<m+1;r+)akr=akr-lki*air; xm-1=am-1m/am-1m-1; for(i=m-2;i>=0;i-)d=0;for( j=i+1;j<m;j+)d=d+aij*xj;xi=(aim-d)/aii;cout<<"该方程组的解为:"<<endl; for(i=0;i<m;i+)cout<<"x"<<i<<"="<<xi<<"t"system("pause");retu

8、rn 0;void load()/*化成三角阵*/*求解*/int i,j;for(i=0;i<m;i+)for( j=0;j<m+1;j+)cin >>aij;运行结果：下面请输入未知数的个数m=3请按顺序输入增广矩阵a:1 2 3 45 1 0 84 6 9 246920-6.5-11.255.50-1.86265e-008-0.1153853.92308该方程组的解为：x0=-9.99999x1=58 x2=-34Press any key to con ti nue总结：列主元素消去法的目的是为了防止减去一个较小的数时大数淹没小数，而使结果产生较大误差，本程序关

9、键在每次消元时找到相应列中的最大项，然后交换两行位置， 在进行计算。(2)LU分解法求解线性方程：#in clude<stdio.h>void solve(float l100,float 叩100,float b,float x,int n) int i,j;float t,s1,s2;float y100;for(i=1;i<=n;i+)/*第一次回代过程开始*/s仁 0;for(j=1;j<i;j+)t=-lij;s1=s1+t*yj;yi=(bi+s1)/lii;for(i=n;i>=1;i-)/*第二次回代过程开始 */s2=0;for( j=n;j&g

10、t;i;j-)t=-ui j;s2=s2+t*xj;xi=(yi+s2)/uii;void mai n()float a100100,l100100,u100100,x100,b100;int i,j,n ,r,k;float s1,s2;1*/for(i=1;i<=99;i+)/*将所有的数组置零，同时将 L矩阵的对角值设为for(j=1;j<=99;j+)lij=O,uij=0;if( j=i) lij=1;printf ("input n:n");/*输入方程组的个数 */sca nf("%d",&n);printf (&quo

11、t;input array A:n");/* 读取原矩阵 A*/for(i=1;i<=n ;i+)for( j=1;j<=n;j+)scan f("%f",&aij);printf ("input array B:n");/*读取列矩阵 B*/for(i=1;i<=n ;i+)scan f("%f", &bi);for(r=1;r<=n;r+)/* 求解矩阵 L 和 U*/ for(i=r;i<=n ;i+)s仁0;for(k=1;k<=r-1;k+)s仁 s1+lrk*u

12、ki; uri=ari-s1;for(i=r+1;i<=n ;i+)s2=0;for(k=1;k<=r-1;k+)s2=s2+lik*ukr;lir=(air-s2)/urr;printf("array L:n");/*输出矩阵 L*/for(i=1;i<=n ;i+)for(j=1;j<=n ;j+)prin tf("%7.3f",lij);prin tf("n");printf("array U:n");/*输出矩阵 U*/for(i=1;i<=n ;i+)for(j=1;j<

13、;=n ;j+)prin tf("%7.3f",uij);prin tf("n");solve(l,u,b,x, n);printf("解为:n");for(i=1;i<=n ;i+)prin tf("x%d=%fn",i,xi);运行结果：in put n:3in put array A:2 2 34 7 7-2 4 5in put array B:3 1 -7array L:1.0000.0000.0002.0001.0000.000-1.0002.0001.000array U:2.0002.0003.

14、0000.0003.0001.0000.0000.0006.000解为：x1=2.000000x2=-2.000000x3=1.000000Press any key to con ti nue总结：关键是把矩阵分解为L、U两个三角矩阵，回代过程比较简单。3:第四章(1)拉格朗日差值多项式；#in clude<stdio.h>#in clude<math.h>#defi ne MAX 100void mai n() int i,j,k,m, n,N,mi;float tmp,mx;float XMAXMAX,YMAX,xMAX,yMAX,aMAX;prin tf(&qu

15、ot;n输入拟合多项式的次数：n");sca nf("%d",&m);printf("n输入给定点的个数 n及坐标(x,y):n");sca nf("%d",&N);prin tf("n");for(i=0;i<N;i+)scan f("%f,%f", &xi, &yi);for(i=0;i<=m;i+)for( j=i;j<=m;j+)tmp=0;for(k=0;k<N;k+)tmp=tmp+pow(xk,(i+j);Xi j=

16、tmp;Xji=Xij;for(i=0;i<=m;i+)tmp=0;for(k=0;k<N;k+)tmp=tmp+yk*pow(xk,i);Yi=tmp;for( j=O;j<m;j+)for(i=j+1,mi=j,mx=fabs(Xjj);i<=m;i+)if(fabs(Xi j)>mx)mi=i;mx=fabs(Xi j);if( j<mi)tmp=Yj;Yj=Ymi;Ymi=tmp;for(k=j;k<=m;k+)tmp=X jk;Xjk=Xmik;Xmik=tmp;for(i=j+1;i<=m;i+)tmp=-Xi j/X jj;Yi+=

17、Yj*tmp;for(k=j;k<=m;k+)Xik+=X jk*tmp;am=Ym/Xmm;for(i=m-1;i>=0;i-)ai=Yi;for( j=i+1;j<=m;j+) ai-=Xij*aj;ai/=Xii;:n");prin tf("n所求的二次多项式为prin tf("P(x)=%f",aO);for(i=1;i<=m;i+)prin tf("+(%f)*xA%d",ai,i);运行结果：输入拟合多项式的次数：5输入给定点的个数 n及坐标（x,y）:31,25,34,2所求的二次多项式为P(x)

18、=1.980417+(0.282759)*xX+(-0.299937)*xA2+(0.022071)*xA3+(0.016624)*xA4+(-0.001934)*xA5Press any key to con ti nue总结：拉格朗日计算公式中，只需要知道各个点即可4:第五章(1)曲线拟合：#in clude<stdio.h>#in clude<math.h>#defi ne MAX 100void mai n()int i,j,k,m, n,N,mi;float tmp,mx;float XMAXMAX,YMAX,xMAX,yMAX,aMAX;prin tf(&q

19、uot;n输入拟合多项式的次数:n");sca nf("%d",&m);printf("n输入给定点的个数 n及坐标(x,y):n");sca nf("%d",&N);prin tf("n");for(i=0;i<N;i+)scan f("%f,%f", &xi, &yi);for(i=0;i<=m;i+)for( j=i;j<=m;j+)tmp=0;for(k=0;k<N;k+)tmp=tmp+pow(xk,(i+j);Xi j

20、=tmp;Xji=Xij;for(i=0;i<=m;i+)tmp=O;for(k=0;k<N;k+)tmp=tmp+yk*pow(xk,i);Yi=tmp;for( j=O;j<m;j+)for(i=j+1,mi=j,mx=fabs(Xjj);i<=m;i+)if(fabs(Xi j)>mx)mi=i;mx=fabs(Xi j);if( j<mi)tmp=Yj;Yj=Ymi;Ymi=tmp;for(k=j;k<=m;k+)tmp=X jk;Xjk=Xmik; Xmik=tmp;for(i=j+1;i<=m;i+)tmp=-Xi j/X jj;Yi

21、+=Yj*tmp; for(k=j;k<=m;k+)Xik+=Xjk*tmp;am=Ym/Xmm;for(i=m-1;i>=0;i-)ai=Yi;for( j=i+1;j<=m;j+) ai-=Xij*aj;ai/=Xii;:n");prin tf("n所求的二次多项式为prin tf("P(x)=%f",aO);for(i=1;i<=m;i+)prin tf("+(%f)*xA%d",ai,i);输入拟合多项式的次数：2输入给定点的个数n及坐标(x,y):51,25.32.48,3-1,5所求的二次多项式为：

22、P(x)=3.952280+(-0.506315)*xA1+(0.050877)*xA2Press any key to continue5:第六章(1)辛普生求积方法：#i nclude <stdio.h>#define N 16/* 等分数 */float fun c(float x) float y;y=4.0/(1+x*x);return(y);void gedia nzhi(float y,float a,float h) int i;for(i=0;i<=N;i+)yi=fu nc(a+i*h);float simps on( float y,float h) f

23、loat s,s1,s2;int i;s1=y1;s2=0.0;for(i=2;i<=N-2;i=i+2) s1+=yi+1;/*计算奇数项的函数值之和*/s2+=yi;/*计算偶数项的函数值之和*/s=y0+yN+4.0*s1+2.0*s2;return(s*h/3.0);main () float a,b,h,s,fN+1;sca nf("%f,%f",&a,&b);h=(b-a)/( float)N;gedia nzhi(f,a,h);s=simps on( f,h);prin tf("s=%fn",s);运行结果：1,3s=

24、1.854590Press any key to con ti nue总结：辛普生算法是一种积分方法，采用三点法插值，如果h较小的话，误差很小，因为4K a它的插值余项 R(f)- - - f()，辛普生算法比较精确，程序关键是对所取的180 2点的取和，注意6:第七章(1)改进欧拉法求解常微分方程的初值问题#i nclude <stdio.h>float fun c(float x,float y) return(y-x);float euler(float xO,float xn,float yO,int N) float x,y,yp,yc,h;int i;x=x0;y=y0

25、;h=(x n-xO)/(float)N;for(i=1;i<=N;i+) yp=y+h*fu nc(x,y);x=x0+i*h;yc=y+h*fu nc(x,yp);y=(yp+yc)/2.0;return(y);main ()float x0,x n, y0,e;int n;prin tf("nin put n :n”)；sca nf("%d",&n);prin tf("i nput xO,x n:n”);sea nf("%f,%f", &xO, &xn);prin tf("i nput y

26、0:n");sea nf("%f", &y0);e=euler(xO,x n,yO,n);prin tf("y(%f)=%6.4f",y0,e); in put n:20in put x0,x n:1,6in put y0:2y(2.000000)=7.0000Press any key to con ti nue(2)四阶龙格一库塔法#i nclude <stdio.h>float fun c(float x,float y) return(x-y);float runge_kutta(float xO,float xn,

27、float yO,int N) float x,y,y1,y2,h,xh;float d1,d2,d3,d4;int i;x=x0;y=y0;h=(x n-xO)/(float)N;for(i=1;i<=N;i+)xh=x+h/2;d1=fu nc(x,y);d2=fu nc(xh,y+h*d1/2.0);d3=fu nc(xh,y+h*d2/2.0);d4=fu nc(xh,y+h*d3);y=y+h*(d1+2*d2+2*d3+d4)/6.0;x=x0+i*h; return(y);main () float x0,x n, y0,e;int N;prin tf("nin put n :n”)；sca nf("%d",&N);prin tf("i nput xO,x n:n”);sea