プログラマ専用SNS ミクプラ

コード:

#include <stdio.h>
#include <math.h>

double func1(double x){
  double y;
  y=x+1/2-sqrt(1-x*x);
  return y;
}

double func2(double x){
  double y;
  y=x+1/2+sqrt(1-x*x);
  return y;
}

int i, j;/*グローバル変数*/

double NewtonRaphson(double (*func_p)(double), double init, double eps, int imax) {/*ニュートン法*/
  double x, x_n, f, df, h;
  x=init;
  i=0;
  h=1e-4;
  while (i<imax) {
    f=func_p(x);
    df=(func_p(x+h)-func_p(x))/h;
    x_n=x - f/df;/*ニュートン反復により次の点x_n を計算する部分*/
    if (fabs(x_n - x)<eps){/* 収束条件を満たしたら解を返す*/
      return x_n;
    }
    i++;/*カウンター*/
    x=x_n;
  }
  return 0.0/0.0; /* NaN を返す*/
}

double Bisection(double (*func_p)(double), double a, double b, double eps){/*二分法*/
  double c;
  c = (a+b)/2;
  j=1;
  while(fabs(a-b)>=2*eps){/*収束前のときループ*/
    c = (a+b)/2;
    j++;/*カウンター*/
    if(func_p(c)>0){/*ｂとf(c)が同符号のとき*/
      c = b;
    }
    if(func_p(c)<0){/*aとf(c)が同符号のとき*/
      c = a;
    }
    if(func_p(c)==0.0){/*f(c)=0のとき*/
      return c;
    }
  }
  return c;
}



int main(void){
  double x;
setbuf(stdout, NULL);
  printf("Newton-Raphson method:\n");/*ニュートン法*/

  printf("Upper part: ");/*関数１*/
  x = NewtonRaphson(func1, 0, 1e-6, 20);
  if (!isnan(x)){
	  printf("(%d times):( %+.6f, %+.6f )\n", i+1, x, func1(x));
  }else{
	  printf("not converged\n");
  }
  x=0;
  printf("Lower part: ");/*関数２*/
  x = NewtonRaphson(func2, 0, 1e-6, 20);
  if (!isnan(x)){
	  printf("(%d times):( %+.6f, %+.6f )\n", i+1, x, func2(x));
  }else{
	  printf("not converged\n");
  }
  x=0;

  printf("Bisection method:\n"); /*二分法*/
  printf("Upper part: "); /*関数１*/

/*ここから下が実行されません*/

  x = Bisection(func1, -1, 1, 1e-6);
  printf("(%d times):( %+.6f, %+.6f )\n", j, x, func1(x));
  x=0.0;
  printf("Lower part: ");/*関数２*/
  x = Bisection(func2, -1, 1, 1e-6);
  printf("(%d times):( %+.6f, %+.6f )\n", j, x, func2(x));
  x=0.0;


return 0;
}

コード:

  y=x-sqrt(1-x*x);
  y=x+sqrt(1-x*x);

コード:

double mysqrt(double target) {
	double left,right,mid,now;
	left=0;right=(target<=1?1:target);
	while(1) {
		mid=(left+right)/2;
		now=mid*mid;
		if(now<target) {
			if(left==mid)break;
			left=mid;
		} else {
			if(right==mid)break;
			right=mid;
		}
	}
	return mid;
}

コード:

double Bisection(double (*func_p)(double), double a, double b, double eps){/*二分法*/
  double c;
  c = (a+b)/2;
  j=1;
  while(fabs(a-b)>=2*eps){/*収束前のときループ*/
    c = (a+b)/2;
    j++;/*カウンター*/
    if(func_p(c)>0){/*ｂとf(c)が同符号のとき*/
      if(c==b)return c;
      c = b;
    }
    if(func_p(c)<0){/*aとf(c)が同符号のとき*/
      if(c==a)return c;
      c = a;
    }
    if(func_p(c)==0.0){/*f(c)=0のとき*/
      return c;
    }
  }
  return c;
}

コード:

    else if(func_p(c)<0){/*aとf(c)が同符号のとき*/
      c = a;
    }

コード:

if(func_p(c)>0){/*ｂとf(c)が同符号のとき*/
  if(c==b)return c;
  c = b;
}
if(func_p(c)<0){/*aとf(c)が同符号のとき*/
  if(c==a)return c;
  c = a;
}

コード:

c = (a+b)/2;