您的位置:首页 > 职场人生

二分法求解方程例子的说明!!

2008-09-13 10:53 393 查看
正在学习C,今天看的例子是二分法求解方程,源代码如下:
#include <stdio.h>
#include <math.h>
#include <malloc.h>
#include <stdlib.h>

double func(double);
int bisectroot(double,double,double,double,double *,int ,int *);

void main()
{
int i,n,m;
double a,b,h,eps,*x;
n=3;
x=(double *)calloc(n,sizeof(double));
if(x==NULL)
{
  printf("memory error!!\n");
  exit(1);
}
a=-1;
b=7;
h=0.1;
eps=1.e-8;
bisectroot(a,b,h,eps,x,n,&n);
printf("the root of y=sin(x) between %2.0f and %2.0f have %d\n",a,b,n);
printf("they are:\n");
for(i=0;i<n;i++)
{
  printf("x[%d]=%e\n",i,x[i]);
}
free(x);
}

double func(double x)
{
return(sin(x));
}

int bisectroot(double a,double b,double h,double eps,double *x,int n,int *m)
{
double z,z0,z1,y,y0,y1;
*m=0;
z=a;
y=func(z);
while(1)
{
  if((z>b+h/2)||(*m==n))
  return(1);
  if(fabs(y)<eps)
  {
  *m+=1;
  x[*m-1]=z;
  z+=h/2;
  y=func(z);
  continue;
  }
  z1=z+h;
  y1=func(z1);
  if(fabs(y1)<eps)
  {
  *m+=1;
  x[*m-1]=z1;
  z=z1+h/2;
  y=func(z);
  continue;
  }
  if(y*y1>0)
  {
  y=y1;
  z=z1;
  continue;
  }
  while(1)
  {
  if(fabs(z1-z)<eps)
  {
  *m+=1;
  x[*m-1]=(z1+z)/2;
  z=z1+h/2;
  y=func(z);
  break;
  }
  z0=(z1+z)/2;
  y0=func(z0);
  if(fabs(y0)<eps)
  {
  *m=*m+1;
  x[*m-1]=z0;
  z=z0+h/2;
  y=func(z);
  break;
  }
  if(y*y0<0)
  {
  z1=z0;
  y1=y0;
  }
  else
  {
  z=z0;
  y=y0;
  }
  }
}
}

可能有的初学者看了这个程序,会以为这个程序使用于任何函数的求解,其实不是这样.
  1.如果有fabs(func(x))和fabs(func(x+h/2))相对于同一点都<eps的情况出现的话,可能会有问题.  2.如果有fabs(func(x))和fabs(func(x+h))在两点的左边且都<eps,再分别加上h/2,形成func(x+h/2)和func(x+h+h/2)在第二个点的两端,就形成了func(x+h/2)*func(x+h+h/2)<0.这样就会重复求得第二个解.   所以,还是要知道要求的函数的大致情况,调整相应的参数.
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签:  职场 休闲 二分法