C语言解题经典程序(2)

totoerer · 发表于 2004-5-14 08:32:21

<

align=center>迭代法求方程的根 DIDAIROT.C
<

>#include<math.h>
<

>#include<stdio.h>
double f(double x)
{ return(2.0-log(x)); }
 
double ca_root(double x0,double eps)
{
int n=0; double x1;
while(n<800)
 {
 x1=f(x0);
 if(fabs(x1-x0)<eps)break;
 x0=x1; n++;
 }
return(x1);
}
main()
{
double root, eps=1.0e-12;
root=ca_root(1.0,eps);
printf("\nThe root =%14.12lf",root);
}经过64次循环，计算结果为：1.557145598997

totoerer · 发表于 2004-5-14 08:36:10

<

center" align=center>二分法求方程的根 ROOT.C<

22pt; mso-line-height-rule: exactly">#include<math.h><

22pt; mso-line-height-rule: exactly">#include<stdio.h>double f(double x){return(x+log(x)-2.0);}double ca_root(double a,double b,double eps){int sg1,sg2,n=1; double x，fa,fx;fa=f(a); sg1=(fa>0.0) ? 1:-1;for(;;) { x=0.5*(a+b); fx=f(x); if(fx==0.0)break; sg2=(fx>0.0) ? 1:-1; if(sg2==sg1)a=x; else b=x; if(fabs(b-a)<eps)break; n++; }printf("\nN=%d ",n); return(x);}main(){double eps=1.0e-12，root;root=ca_root(1.0,2.0,eps);printf("\nThe Root=%14.12lf",root);}经过40次循环，计算结果为：1.557145598998

totoerer · 发表于 2004-5-14 08:37:30

<

center" align=center>斐波那契法求函数的极大值fibona.c<

22pt; mso-line-height-rule: exactly">#include<math.h><

22pt; mso-line-height-rule: exactly">#include<stdio.h>#include<stdlib.h> double f(double x){return(x*(2.0-log(x)));} double fibona(double a,double b,double eps){int n=1,i;double c,h,x1,x2,f1=1.0,f2,f0=1.0;c=(b-a)/eps;for(;;) { f2=f0+f1; if(f2>=c)break; f0=f1;f1=f2; n++; }c=(b-a)/f2;h=f0*c;x1=a+h;f1=f(x1);x2=b-h;f2=f(x2);for(i=n;i>=2;i--) { if(f1>=f2) { b=x2;h=x2-x1;x2=x1; f2=f1;x1=a+h;f1=f(x1); } else { a=x1;h=x2-x1;x1=x2;f1=f2; x2=b-h;f2=f(x2); } if(fabs(x2-x1)<eps)break; }return((a+b)*0.5);} main(){double mx=0.0,my=0.0;mx=fibona(1.0,4.0,1.0e-8);my=f(mx);printf("\n Maxx=%14.12lf Maxy=%14.12lf",mx,my);getch();}计算结果：Maxx=2.718281828764Maxy=2.718281828459

totoerer · 发表于 2004-5-14 08:43:10

<

center" align=center>变步长改进欧拉法求微分方程的数值解 EULER.C<

24pt; mso-line-height-rule: exactly">#include<math.h><

24pt; mso-line-height-rule: exactly">#include<stdio.h> double f(double x,double y){ return(x*exp(-y)); } double euler(double x0,double y0,double xn,unsigned long n){unsigned long k; double t,x,y,y1,h;h=(xn-x0)/(double)n; x=x0; y=y0;for(k=1;k<=n;k++) { t=f(x,y); y1=y+h*t; y1=y+h*0.5*(t+f(x+h,y1)); x+=h;y=y1; }return(y1);} double buzh(double x0,double y0,double xn,double eps){int k=0; unsigned long m=64;double y,y1;y=euler(x0,y0,xn,m);while(k<20) { m*=2; y1=euler(x0,y0,xn,m); if(fabs(y1-y)<eps)break; y=y1;k++; printf("\n k=%u y=%15.14lf",k,y1); }return(y1);} main(){double x0=0.0,y0=0.0,answer,eps=1.0e-012;printf("\nThe Answer is %15.14lf",log(13.5));answer=buzh(x0,y0,5.0,eps);printf("\nThe Answer is %15.14lf",answer);} 循环14次，计算结果：<v:shapetype> <v:stroke joinstyle="miter"></v:stroke><v:formulas><v:f eqn="if lineDrawn pixelLineWidth 0"></v:f><v:f eqn="sum @0 1 0"></v:f><v:f eqn="sum 0 0 @1"></v:f><v:f eqn="prod @2 1 2"></v:f><v:f eqn="prod @3 21600 pixelWidth"></v:f><v:f eqn="prod @3 21600 pixelHeight"></v:f><v:f eqn="sum @0 0 1"></v:f><v:f eqn="prod @6 1 2"></v:f><v:f eqn="prod @7 21600 pixelWidth"></v:f><v:f eqn="sum @8 21600 0"></v:f><v:f eqn="prod @7 21600 pixelHeight"></v:f><v:f eqn="sum @10 21600 0"></v:f></v:formulas><v:path connecttype="rect" gradientshapeok="t" extrusionok="f"></v:path><lock aspectratio="t" v:ext="edit"></lock></v:shapetype><v:shape><v:imagedata></v:imagedata></v:shape>

totoerer · 发表于 2004-5-14 08:44:54

<

>/*改进坐标轮换法求多元函数极小值 DIRECT.C */<

>#include<conio.h><

>#include<math.h>#include<stdio.h>#define N 2int j; double x[N];double fun(double x[]){return(cos(x[0]*x[0])+x[0]-x[1]+2.0*x[0]*x[0]+2.0*x[0]*x[1]+x[1]*x[1]+log(1.0+x[1]*x[1]));}double fx(double y[],double lbt){int i;for(i=0;i<N;i++) x=y; x[j]=y[j]+lbt;return(fun(x));} double fy(double y[],double lbt){int i; double x2[2];for(i=0;i<N;i++) x2=x+lbt*y;return(fun(x2));}double ca_min(double a,double b,double eps,double (*fun)(),double y[]){double c,x1,x2,x3,f1,f2,f3;c=(b-a)*0.25; x1=a+c; f1=fun(y,x1);x2=x1+c; f2=fun(y,x2); x3=x2+c; f3=fun(y,x3);for(;;) { if(f1<=f2&&f1<f3) { b=x2; x2=x1; f2=f1; } else if(f3<=f2&&f3<f1) { a=x2; x2=x3; f2=f3; } else { a=x1; b=x3; } c=(b-a)*0.25; if(fabs(c)<eps)break; x1=a+c; f1=fun(y,x1); x3=x2+c; f3=fun(y,x3); }return(x2);}int direct(double x0[],double eps){int i; unsigned k=0; double lbt,px,x1[N],y[N];for(i=0;i<N;i++)y=x0;while(k<200) { for(j=0;j<N;j++) { lbt=ca_min(-1.0,1.0,eps*0.01,fx,y); y[j]=y[j]+lbt; } px=0.0; for(i=0;i<N;i++) { x=x1=y; y=x1-x0; x0=x1; px+=(y*y); } if(px<eps)break; lbt=ca_min(-2.0,2.0,eps*0.1,fy,y); for(i=0;i<N;i++) y=x1+lbt*y; k++; }return(1);} main(){int i;double x0[N],eps=1.0e-014;x0[0]=2.0;x0[1]=3.0;direct(x0,eps);printf("\nf(X) = %12.10lf ",fun(x0)); for(i=0;i<N;i++) printf("\nx[%d] = %12.10lf ",i,x0);}

totoerer · 发表于 2004-5-14 08:45:33

<

18pt; TEXT-ALIGN: center; mso-line-height-rule: exactly" align=center>黄金分割法求函数的极大值gold1.c<

18pt; mso-line-height-rule: exactly">#include<math.h><

18pt; mso-line-height-rule: exactly">#include<stdio.h>double f(double x){return(x*sin(x));}double ca_max(double a,double b,double eps){int n=0;double c,x1,x2,f1,f2;c=0.5*(sqrt(5.0)-1);x1=a+c*(b-a); f1=f(x1); x2=a+b-x1; f2=f(x2);for(;;) { if(f1<f2) { b=x1; x1=x2; f1=f2; x2=b-c*(b-a); f2=f(x2); } else if(f1>f2) { a=x2; x2=x1; f2=f1; x1=a+c*(b-a); f1=f(x1); } else { a=x2; b=x1; x1=a+c*(b-a); x2=a+b-x1; f1=f(x1); f2=f(x2); } if(fabs(x1-x2)<eps)break; n++; }c=(x1+x2)*0.5; printf("\nn=%d",n);return(c);}main( ){double mx,my;mx=ca_max(0.0,3.0,1.0e-10); my=f(mx);printf("\n Maxx=%14.12lf Maxy=%14.12lf",mx,my);}计算结果：循环次数：n=43极大值点：Maxx=2.028757834797极大值： Maxy=1.819705741160

totoerer · 发表于 2004-5-14 08:46:24

<

1.2pt 0cm">程序一、计算圆周率π<

>#include<stdio.h><

>main( ) { double a,m=3.0,pi,n=1.0; a=1.0/3.0;pi=1.0+a; do { n+=1.0; m+=2.0; a=a*n/m; pi+=a; }while(a>1.0e-10); printf("\npi=%14.12lf ",pi*2.0); }

totoerer · 发表于 2004-5-14 08:46:46

<

1.2pt 0cm">程序二、计算任意正数的平方根<

>#include<stdio.h><

>#include<math.h>main( ) { double a=3.0,x0=1.0,x1,b,eps=1.0e-10; printf(“Input a=?”); scanf(“%lf”,&a); for(;;) { x1=(x0+a/x0)*0.5; if(fabs(x1-x0)<eps)break; x0=x1; } printf("\na=%g The answer is %14.12lf",x1); }

totoerer · 发表于 2004-5-14 08:53:23

<

center" align=center>简化牛顿法解非线性方程组newton1.c<

23pt; mso-line-height-rule: exactly">#include<math.h><

23pt; mso-line-height-rule: exactly">#include<stdio.h>#define N 3 double pi,a[N][N],b[N][N],fx[N]; void F(double x[]){double ex,sx;fx[0]=3.0*x[0]-cos(x[1]*x[2])-0.5;fx[1]=x[0]*x[0]-81.0*(x[1]+0.1)*(x[1]+0.1)+sin(x[2])+1.06;ex=exp(-x[0]*x[1]);fx[2]=ex+20.0*x[2]+pi*10.0/3.0-1.0;} void DF(double x[]){double ex,sx;ex=exp(-x[0]*x[1]);a[0][0]=3.0;sx=sin(x[1]*x[2]);a[0][1]=x[2]*sx;a[0][2]=x[1]*sx;a[1][0]=2.0*x[0];a[1][1]=-162.0*(x[1]+0.1);a[1][2]=cos(x[2]);a[2][0]=-x[1]*ex;a[2][1]=-x[0]*ex;a[2][2]=20.0;}int inverta(int n) { int i,j,k,i0; double c,t,maxa,eps=1.0e-15; for(i=0;i<n;i++) for(j=0;j<n;j++) if(j==i)b[j]=1.0; else b[j]=0.0; for(k=0;k<n;k++) { maxa=0.0; for(i=k;i<n;i++) { t=fabs(a[k]); if(t>maxa) { maxa=t;i0=i; } } if(maxa<=eps)return(0); if(i0!=k) { for(j=0;j<n;j++) { t=a[k][j]; c=b[k][j]; a[k][j]=a[i0][j]; b[k][j]=b[i0][j]; a[i0][j]=t; b[i0][j]=c; } } c=1.0/a[k][k]; for(j=k;j<n;j++)a[k][j]*=c; for(j=0;j<n;j++)b[k][j]*=c; for(i=k+1;i<n;i++) { for(j=k+1;j<n;j++) a[j]-=(a[k]*a[k][j]); for(j=0;j<n;j++) b[j]-=(a[k]*b[k][j]); } } for(k=n-2;k>=0;k--) for(i=k;i>=0;i--) for(j=0;j<n;j++)b[j]-=(b[k+1][j]*a[k+1]); return(1); } int newtonfa(double x0[],double eps){int i,j,k=1;double mx;double x1[N],dx[N],y[N];DF(x0);if(!inverta(N))return(0);for(;;) { F(x0); for(i=0;i<N;i++) { dx=0.0; for(j=0;j<N;j++) dx+=b[j]*fx[j]; } mx=0.0; for(i=0;i<N;i++)if(fabs(dx)>mx)mx=fabs(dx); for(i=0;i<N;i++)x1=x0-dx; for(i=0;i<N;i++)x0=x1; if(mx<eps)break; printf("\nk=%d",k); printf("\nx1=%12.10lf x2=%12.10lf x3=%12.10lf",x1[0],x1[1],x1[2]); k++;if(k>400)return(0); }return(1);} main(){int i;double x0[N]={2.0,1.0,0.0},eps=1.0e-14;pi=4.0*atan(1.0);if(newtonfa(x0,eps)) { for(i=0;i<N;i++) printf("\nx[%d] = %14.12lf ",i,x0); }else printf("\nFiled !");}循环293次，计算结果：<v:shapetype> <v:stroke joinstyle="miter"></v:stroke><v:formulas><v:f eqn="if lineDrawn pixelLineWidth 0"></v:f><v:f eqn="sum @0 1 0"></v:f><v:f eqn="sum 0 0 @1"></v:f><v:f eqn="prod @2 1 2"></v:f><v:f eqn="prod @3 21600 pixelWidth"></v:f><v:f eqn="prod @3 21600 pixelHeight"></v:f><v:f eqn="sum @0 0 1"></v:f><v:f eqn="prod @6 1 2"></v:f><v:f eqn="prod @7 21600 pixelWidth"></v:f><v:f eqn="sum @8 21600 0"></v:f><v:f eqn="prod @7 21600 pixelHeight"></v:f><v:f eqn="sum @10 21600 0"></v:f></v:formulas><v:path gradientshapeok="t" connecttype="rect" extrusionok="f"></v:path><lock v:ext="edit" aspectratio="t"></lock></v:shapetype><v:shape><v:imagedata></v:imagedata></v:shape>，<v:shape> <v:imagedata></v:imagedata></v:shape>，<v:shape> <v:imagedata></v:imagedata></v:shape>

totoerer · 发表于 2004-5-14 08:53:26

<

center" align=center>简化牛顿法解非线性方程组newton1.c<

23pt; mso-line-height-rule: exactly">#include<math.h><

23pt; mso-line-height-rule: exactly">#include<stdio.h>#define N 3 double pi,a[N][N],b[N][N],fx[N]; void F(double x[]){double ex,sx;fx[0]=3.0*x[0]-cos(x[1]*x[2])-0.5;fx[1]=x[0]*x[0]-81.0*(x[1]+0.1)*(x[1]+0.1)+sin(x[2])+1.06;ex=exp(-x[0]*x[1]);fx[2]=ex+20.0*x[2]+pi*10.0/3.0-1.0;} void DF(double x[]){double ex,sx;ex=exp(-x[0]*x[1]);a[0][0]=3.0;sx=sin(x[1]*x[2]);a[0][1]=x[2]*sx;a[0][2]=x[1]*sx;a[1][0]=2.0*x[0];a[1][1]=-162.0*(x[1]+0.1);a[1][2]=cos(x[2]);a[2][0]=-x[1]*ex;a[2][1]=-x[0]*ex;a[2][2]=20.0;}int inverta(int n) { int i,j,k,i0; double c,t,maxa,eps=1.0e-15; for(i=0;i<n;i++) for(j=0;j<n;j++) if(j==i)b[j]=1.0; else b[j]=0.0; for(k=0;k<n;k++) { maxa=0.0; for(i=k;i<n;i++) { t=fabs(a[k]); if(t>maxa) { maxa=t;i0=i; } } if(maxa<=eps)return(0); if(i0!=k) { for(j=0;j<n;j++) { t=a[k][j]; c=b[k][j]; a[k][j]=a[i0][j]; b[k][j]=b[i0][j]; a[i0][j]=t; b[i0][j]=c; } } c=1.0/a[k][k]; for(j=k;j<n;j++)a[k][j]*=c; for(j=0;j<n;j++)b[k][j]*=c; for(i=k+1;i<n;i++) { for(j=k+1;j<n;j++) a[j]-=(a[k]*a[k][j]); for(j=0;j<n;j++) b[j]-=(a[k]*b[k][j]); } } for(k=n-2;k>=0;k--) for(i=k;i>=0;i--) for(j=0;j<n;j++)b[j]-=(b[k+1][j]*a[k+1]); return(1); } int newtonfa(double x0[],double eps){int i,j,k=1;double mx;double x1[N],dx[N],y[N];DF(x0);if(!inverta(N))return(0);for(;;) { F(x0); for(i=0;i<N;i++) { dx=0.0; for(j=0;j<N;j++) dx+=b[j]*fx[j]; } mx=0.0; for(i=0;i<N;i++)if(fabs(dx)>mx)mx=fabs(dx); for(i=0;i<N;i++)x1=x0-dx; for(i=0;i<N;i++)x0=x1; if(mx<eps)break; printf("\nk=%d",k); printf("\nx1=%12.10lf x2=%12.10lf x3=%12.10lf",x1[0],x1[1],x1[2]); k++;if(k>400)return(0); }return(1);} main(){int i;double x0[N]={2.0,1.0,0.0},eps=1.0e-14;pi=4.0*atan(1.0);if(newtonfa(x0,eps)) { for(i=0;i<N;i++) printf("\nx[%d] = %14.12lf ",i,x0); }else printf("\nFiled !");}循环293次，计算结果：<v:shapetype> <v:stroke joinstyle="miter"></v:stroke><v:formulas><v:f eqn="if lineDrawn pixelLineWidth 0"></v:f><v:f eqn="sum @0 1 0"></v:f><v:f eqn="sum 0 0 @1"></v:f><v:f eqn="prod @2 1 2"></v:f><v:f eqn="prod @3 21600 pixelWidth"></v:f><v:f eqn="prod @3 21600 pixelHeight"></v:f><v:f eqn="sum @0 0 1"></v:f><v:f eqn="prod @6 1 2"></v:f><v:f eqn="prod @7 21600 pixelWidth"></v:f><v:f eqn="sum @8 21600 0"></v:f><v:f eqn="prod @7 21600 pixelHeight"></v:f><v:f eqn="sum @10 21600 0"></v:f></v:formulas><v:path gradientshapeok="t" connecttype="rect" extrusionok="f"></v:path><lock v:ext="edit" aspectratio="t"></lock></v:shapetype><v:shape><v:imagedata></v:imagedata></v:shape>，<v:shape> <v:imagedata></v:imagedata></v:shape>，<v:shape> <v:imagedata></v:imagedata></v:shape>

		自动登录	找回密码
密码			注-册-帐-号