fork(4) download
copy 
  1. #include <stdio.h>
  2. #include <math.h>
  3. #include <iostream>
  4. #include<cfloat>
  5. using namespace std;
  6.  
  7. const int max_iterations=100;
  8. const double interval = 0.001;  // interval length for roots subdivision
  9. bool is_root = false;  // is false if function has a second type break in root candidate point
  10. double delta; //precision of root search
  11. double A, B;
  12.  
  13. double f(double x)//our function
  14. {
  15.  double r = DBL_MAX;
  16.  try
  17.   {
  18.    r =  x / (2*sin(x)+1) - tan(log(x*x+1));  
  19.   }
  20.  catch(...)
  21.   {
  22.    //nothing to do for arithmetic exception
  23.   }
  24.  return r;
  25. }
  26.  
  27. bool is_change_sign(double fa, double fb)
  28. {
  29.  return fa*fb<0;
  30. }
  31.  
  32. double find_root(double a, double b, double delta,  bool *yes)
  33. {
  34.  double fa;
  35.  double fb;
  36.  double fc;
  37.  double c;
  38.  int l=0;
  39.  fa=f(a);
  40.  fb=f(b);
  41.  double dy,dy1; //old and new OY projections of chord
  42.  dy=fabs(fa-fb);
  43.  dy1=DBL_MAX;
  44.  while(true)
  45.    {
  46. 	l++;
  47. 	if(l > max_iterations)
  48. 	  {
  49. 	   *yes=false;
  50. 	   return (a+b)/2;
  51. 	  }
  52. 	c = a - (b - a) * fa/(fb - fa);
  53. 	fc = f(c);
  54. 	if(fa*fc>0)
  55. 	  {
  56. 	   a=c;
  57. 	   fa=fc;
  58. 	  }
  59. 	else
  60. 	  {
  61. 	   b=c;
  62. 	   fb=fc;
  63. 	  }
  64. 	if(fabs(fa)<delta)
  65. 	  {
  66. 	   *yes=true;
  67. 	   return a;
  68. 	  }
  69. 	if(fabs(fb)<delta)
  70. 	  {
  71. 	   *yes=true;
  72. 	   return b;
  73. 	  }
  74. 	if(l>10) //there's a magic number, we must check conditions of break or stop after some iterations(some = 10)
  75. 	  {
  76. 	   if(dy1>dy) //checking if there's break of second kind(pole)
  77. 	     {
  78. 	      *yes=false;
  79. 	      return (a+b)/2;
  80. 	     }
  81. 	   if(fabs(b-a)<delta) //else we are near the root
  82. 		 {
  83. 		  *yes=true;
  84. 		  return (a+b)/2;
  85. 		 }
  86. 	  }
  87. 	dy=dy1; //change old chord height
  88. 	dy1=fabs(fa-fb); //calculate new chord height
  89.    }
  90.  // we can't get here
  91.  *yes=false;
  92.  return c;
  93. }
  94.  
  95. int main()
  96. {
  97.  scanf("%lf %lf %lf",&A,&B,&delta);
  98.  double a,b;
  99.  double fa, fb;
  100.  /*
  101.  test interval for finding roots
  102.  int k=2;
  103.  A= sqrt(exp(M_PI*(1/2.0+k))-1)-M_PI;//70010637.8772-10;
  104.  B= sqrt(exp(M_PI*(1/2.0+k))-1)+M_PI;//70010637.8772+10;
  105.  delta = 1e-12;
  106.  printf("A=%20.19g A=%20.19g delta=%20.19g \n",A,B,delta);
  107.  */
  108.  fb = f(A);
  109.  for(a = A, b = A + interval; b<=B; a += interval, b += interval)
  110.   {
  111.    fa = fb;
  112.    fb = f(b);
  113.    if(is_change_sign(fa,fb))
  114.      {
  115.       double r=find_root(a,b,delta, &is_root);
  116.       if(is_root)
  117.        {
  118.         printf("%14.15lf \n",r);
  119.        }
  120.      }
  121.   }
  122.  return 0;
  123. }
  124.  
Success #stdin #stdout 0s 3460KB
comments (?)
stdin
copy 
-10 10 0.0001
stdout
copy 
-1.257414829166168 
-0.000000000000103 
0.547316038071892
https://ideone.com/uONVYq
Mloop2
language:
C++14 (clang 8.0)
created:
9 years ago
visibility:
secret

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!