数值分析计算方法程序汇总.doc

《数值分析计算方法程序汇总.doc》由会员分享，可在线阅读，更多相关《数值分析计算方法程序汇总.doc（13页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、【精品文档】如有侵权，请联系网站删除，仅供学习与交流数值分析计算方法程序汇总.精品文档.（一）秦九韶法例：已知p=2x2-3x-4，求x=1时，p=?程序：#include stdio.hvoid main()float a50,b,x;int n,i,k;scanf(%d%f,&n,&x);for(i=0;i=n;i+)scanf(%f,&ai);b=a0;k=1;while(k=n)b=x*b+ak;k=k+1;printf(%f,b);结果：-5.000000(二）复化中矩形公式例：求=g(y1)+g(y2)+g(ym)h的值(已知g(y)=yn/2-1e-y/2 ，n=11，h=0.1

2、，ym=20)。程序：#include stdio.h#include main() double y,h,s;int n,i;n=11;h=0.1;s=0;y=0;for(i=0;i=200;i+) y=y+h; s=s+pow(y,n/2.0-1)*exp(-y/2.0);s=s*h;printf(%fn,s);结果：2266.139467（三）复化梯形公式例：求程序：#include stdio.hvoid main()double a,b,s,h,x;int i,n;a=-1.0;b=1.0;n=10;h=(b-a)/n;x=a;s=x*x/2;for(i=1;in;i+)x=x+h;

3、s=s+x*x;s=s+b*b/2;s=s*h;printf(s=%fn,s);结果：s=0.680000（四）复化辛普森公式例：求程序：#include stdio.hvoid main()double a,b,c,s,h,x,y;int i,n;a=0.0;b=1.0;n=10;s=0.0;h=(b-a)/n;x=a;y=x+h;c=(x+y)/2;for(i=1;i=n;i+)s=s+x*x*x+4*c*c*c+y*y*y;x=x+h;y=y+h;c=c+h;s=s*h/6;printf(s=%fn,s);结果：s=0.250000（五）复化高斯公式例：求程序：#include #inc

4、lude main()double a,b,h,s,x1,x2; int i,n; a=0;b=2;n=20;s=0; h=(b-a)/n; for(i=0;in;i+) x1=a+i*h+h/2*(1/1.732+1); x2=a+i*h+h/2*(1-1/1.732); s=s+x1*x1*x1+x2*x2*x2; s=h/2*s; printf(s=%fn,s);结果：s=4.000000(六）二维中矩形公式例：求程序：#include #include main()double a,b,c,d,x50,y50,hx,hy,s; int i,j,n,k; a=0.0;b=1.0;c=0.

5、0;d=1.0; n=k=10; hx=(b-a)/n;hy=(d-c)/k; x0=a+hx/2;y0=c+hy/2;s=0.0; for(j=0;jk;j+) for(i=0;in;i+) s=s+xi*xi+yj*yj; xi+1=xi+hx; yj+1=yj+hy; s=hx*hy*s; printf(s=%fn,s);结果：s=0.665000(七）迭代法例：求x=x2的解。程序：#include stdio.h#includemain() double x,xl,y,yl;int i,j;x=0.5;xl=x;y=0.5;yl=y;for(i=0;i+)x=x*x; if(fabs

6、(xl-x)0.0001) break; else xl=x;for(j=0;j+)y=sqrt(y); if(fabs(yl-y)0.0001) break; else yl=y;printf(x=%f,y=%fn,x,y);结果：x=0.000000，y=0.999915（八）牛顿迭代法y=f(x),求f(x*)=0。y-y0=f(x0)(x-x0)-y0=f(x0)(x1-x0)x1=x0-f(x0)/f(x0)xn+1=xn-f(xn)/f(xn)例：求f(x)=x-x2的零点(已知x0=0.1)。程序：#include stdio.h#includemain() double x,x

7、l;int i;x=0.1;xl=x;for(i=0;i+)x=x-(x-x*x)/(1-2*x); if(fabs(xl-x)0.0001) break; else xl=x;printf(x=%f,x);结果：x=-0.000000（九）二分法已知f(a)f(b)0，则在a,b中至少有一点f(x*)=0。令x1=(a+b)/2如f(a)f(x1)0,取a1=a,b1=x1；否则，取a1=x1,b1=b。令x2=(a1+b1)/2如f(a1)f(x2)0,取a2=a1,b2=x2；否则，取a2=x2,b2=b1。令xn+1=(an+bn)/2如f(an)f(xn+1)0,取an+1=an,b

8、n+1=xn+1；否则，取an+1=xn+1,bn+1=bn。例：求f(x)=x-x2在0.1,2.45内的零点。程序：#include stdio.h#includemain() double a,al,b,bl,x,xl; int i;a=0.1;b=2.45;for(i=0;i+)x=(a+b)/2; if(a-a*a)*(x-x*x)0) al=a,b=x; else a=x,bl=b; if(fabs(xl-x)0.0001) break; else xl=x;printf(x=%f,xl);结果：x=1.000040（十）一元回归已知（x1,y1),(x2,y2)(xN,yN),y

9、=a+bx yi=a+bxi+Ei Q(a,b)= Ei=yi-a-bxi总误差Q= 令得aN+bX-Y=0, aX+bX2-XY=0（其中X=, Y=, X2=, XY=）。可得：例：已知（0，0），（1，2），（2，5）；求a，b。程序：#include stdio.hmain() float a,b,E,X,Y,Xl,Xy,Q,x3=0.0,1.0,2.0,y3=0.0,2.0,5.0; int i,j,n=3; X=0;Y=0;Xl=0;Xy=0; for(i=0;in;i+) X=X+xi; Y=Y+yi; Xl=Xl+xi*xi; Xy=Xy+xi*yi; a=(Xl*Y-X*Xy

10、)/(n*Xl-X*X); b=(n*Xy-X*Y)/(n*Xl-X*X); printf(a=%f,b=%fn,a,b); Q=0;E=0; for(j=0;j=3;j+) Q=Q+E*E; E=yj-a-b*xj; printf(Q(总误差)=%fn,Q);结果：a=-0.166667,b=2.500000Q(总误差)=0.166667（十一）泰勒插值例：计算当x=2.104时，ex=pn(x)=？程序：#include stdio.h#include void main()double x,l,p; int i,n; x=2.104;n=100; p=1;l=1; for(i=0;i=n

11、;i+) l=l*x/(i+1); p=p+l; printf(p=%fn,p);结果：p=8.198900（十二）拉格朗日插值例：已知（0，1），（3，2），（8，3），计算x=5时，求p1(x)=？程序：#include stdio.h#include void main()float x3=0.0,3.0,8.0,y3=1.0,2.0,3.0,xl,p,l; int n,k,j; n=2; xl=5.0;p=0.0; for(k=0;k=n;k+) l=1; for(j=0;j=n;j+) if(j=k) continue; else l=l*(xl-xj)/(xk-xj); p=p+l

12、*yk; printf(p=%fn,p);结果：p=2.500000（十三）牛顿插值已知（x0,y0),(x1,y1),(xn,yn);求pn=a0+a1x+a2x2+anxn（1）将上式依次代入（1）中例：已知（0，0），（1，1），（2，4）；求p2(x)=？程序：#include #include main()double x3=0,1,2,y3=0,1,4,g3,gl,p,pl,X; int i,n=2;X=4; g0=y0; g1=(y1-y0)/(x1-x0); gl=(y2-y0)/(x2-x0); g2=(gl-g1)/(x2-x1); pl=gn; for(i=n;i=1;i

13、-) p=pl*(X-xi-1)+gi-1; pl=p; printf(p=%fn,p);结果：p=16.000000（十四）分段插值例：程序：#include stdio.h#include void main()double x10,y10,p10,h,xl;int i,n;n=10;x0=0.0;x9=2.0;xl=2.0;h=(x9-x0)/n;scanf(%d,&i);xi=h*i;yi=exp(xi);xi+1=h*(i+1);yi+1=exp(xi+1);pi+1=yi+(yi+1-yi)/(xi+1-xi)*(xl-xi);printf(pi+1=%fn,pi+1);结果：输入

14、i=0时，pi+1=3.214028（十五）解线性方程组例：解方程组程序：#include stdio.hvoid main()double a33=5.0,1.0,1.0,1.0,4.0,-1.0,1.0,-1.0,6.0,b3=1.0,2.0,3.0,x3; int i,j,k,n; n=3; for(k=0;kn-1;k+) for(j=k+1;jn;j+) for(i=k+1;i=0;k-) xk=(bk-akk*xk+1-akn*xn)/akk; printf(%fn,xk); printf(%fn,xn-1);结果：-0.414921 0.414921 0.572816（十六）列主

15、元消去法例：解方程组程序：#include#include #define N 100 float aNN+1;void main( ) int i,j,k,n; float t,s=0; printf(输入矩阵阶数:); scanf(%d,&n); printf(n); printf(输入增广矩阵:n); for(i=0;in;i+) for(j=0;jn+1;j+) scanf(%f,&aij); for(k=0;kn-1;k+) for(i=k+1;i fabs(akk) ) for(j=k;jn+1;j+) t=akj; akj=aij; aij=t; for(i=k+1;in;i+)

16、 aik=aik / akk; for(j=k+1;j=0;k-) s=0; for(j=k+1;jn;j+) s+=akj*ajn; akn=( akn-s ) / akk; for(i=0;in;i+) printf( x%d=%fn,i+1,ain);结果：x1=2.000000 x2=1.000000 x3=-0.000000（十七）追赶法例：程序：#include stdio.hvoid main()double x3,a4=0,1,1,1,b4=4,4,4,4,c4=1,1,1,0,d4=1,2,1,0,e4,f4; int k,n=3; e0=d0/b0;f0=c0/b0; fo

17、r(k=1;k=0;k-) xk=ek-fk*xk+1; printf(%fn,xk);结果：-0.036900 0.147601 0.377035 0.155741（十八）高斯-赛德尔迭代法例：程序：#include stdio.h#include #define N 100void main()double aNN=8,-3,2,4,11,-1,6,3,12,bN=20,33,36,xN=0,0,0,xl,s; int i,j,k,n=3; for(k=1;k=10;k+) for(i=0;in;i+) s=0.0; xl=xi; for(j=0;jn;j+) if(j=i) contin

18、ue; else s=s+aij*xj; xi=(bi-s)/aii; for(i=0;in;i+) printf( x%d=%fn,i+1,xi);结果：x1=3.000000x2=2.000000x3=1.000000（十九)一步欧拉公式例：已知求x=1时的值。程序：#include stdio.hvoid main()double x,y,h; int i,n=10;h=0.1;x=0;y=1;for(i=1;i=n;i+)y=y+h*(x+y); x=x+h;printf(%fn,y);结果：3.187485（二十）隐式欧拉公式：例：已知求x=1时的值。程序：#include stdi

19、o.h#include void main()double x10,y10,yl,h,e; int i,j,n=10;h=0.1;e=0.01;x0=0.0;y0=1.0;for(i=0;in;i+) yi+1=yi+h*(xi+yi); xi+1=xi+h; for(j=0;j+) yj+1=yj+h*(xj+1+yj+1); yl=yj+1; if(fabs(yl-yi+1)e) break; else yj+1=yl;printf(%fn,yl);结果：3.728783（二十一）四阶龙格-库塔公式例：已知y=x+y，y(0)=1；求y(0.2)=？程序：#include stdio.hv

20、oid main()double x,y,h,k1,k2,k3,k4;int i,n;h=0.1;n=2;x=0.0;y=1.0;for(i=0;in;i+)k1=x+y;k2=x+y+h/2+h*k1/2;k3=x+y+h/2+h*k2/2;k4=x+y+h+h*k3;y=y+h*(k1+2*k2+2*k3+k4)/6;x=x+h;printf(y(0.2)=%fn,y);结果：y(0.2)=1.242805补充例题：（二十二）例：已知y=x+y+y，用一步欧拉公式求n=10时，y=？，y=？程序：#include stdio.hvoid main()double x,y,z,h;int i

21、,n;h=0.1;n=10;x=0;y=0;z=1;for(i=1;i=n;i+)y=y+h*z;z=z+h*(x+y+z);x=x+h;printf(y=%f,y=%fn,y,z);结果：y=2.001568,y=4.314725（二十三）例：已知y=x+y，y(0)=1，y(1)=0，求n=3时，y(1)=？y(2)=？程序：#include stdio.hvoid main()double x3,y3,d3,e3,f3,b,h;int i,k,n=3;h=1.0/n;b=-2.0-h*h;x1=h;d1=x1*h*h-1.0;for(i=2;in;i+)xi=xi-1+h;di=h*h*xi;e1=d1/b;f1=-1/b;for(k=2;k=1;k-) yk=ek+fk*yk+1; printf(y%d=%f,k,yk);结果：y2=0.233333y1=0.566667