Unable to increase the precision in computing square root using bisection method.

[ajiten73]

Tried to code bisection method to compute square root of 2, till the difference between two consecutive steps is less than some specific value, say 10^{-11}, where tolerance variable is specified as 11.

But, inspite of choosing long double data type,am unable to get beyond 6 places after decimal.

The same result is obtained, if choose tolerance as 11, 30, 5.

Cannot understand the reason for that.

The reason for that is unclear, as know that can get infinite number of decimal digits in the decimal expansion of sqrt{2}.

#include <iostream>
#include <math.h>
using namespace std;

int main() {
   int val, tolerance;
   cout << "Input the values of the number whose square root need be found, and the tolerance limit in the number of digits, of integer variable: tolerance, giving tolerance limit as: 10^{-tolerance}" << endl;
   cout<< "Enter the number, followed by the # tolerance digits."<< endl;
   cin >>val>> tolerance;
   if(val< 0)
      exit(0); 
   else 
   {
      long double low = 0, high = val;
      long double  mid = (low + high)/2;
      int c = 0; 
      long double diff = 0;

      while (!c) 
      {
         if((mid*mid-val)<0)
             diff = val - mid*mid;
         if (diff <= pow(10, -1*tolerance))
         {
            c = 1;
            cout << "Square root of 2= "<< mid<<endl;
         }
         else if(mid * mid > val) 
         {
               high = mid;
               mid = (low + high)/2;      
          } 
          else 
          {
               low = mid;
               mid = (low + high)/2;
          }
      }
    }
   return 0;
}

Am running on onlinegdb.com, using the C++ option.

[salem c]

A couple of changes.

#include <iostream>
#include <iomanip>
#include <math.h>
using namespace std;

int main() {
   int val, tolerance;
   cout << "Input the values of the number whose square root need be found, and the tolerance limit in the number of digits, of integer variable: tolerance, giving tolerance limit as: 10^{-tolerance}" << endl;
   cout<< "Enter the number, followed by the # tolerance digits."<< endl;
   cin >>val>> tolerance;
   if(val< 0)
      exit(0); 
   else 
   {
      long double limit = powl(10,-tolerance);  // use powl, not pow for better precision
      long double low = 0, high = val;
      long double  mid = (low + high)/2;
      int c = 0; 
      long double diff = 0;

      while (!c) 
      {
         if((mid*mid-val)<0)
             diff = val - mid*mid;
         if (diff <= limit)
         {
            c = 1;
            cout << "Square root of 2= "<< setprecision(tolerance) << mid<<endl;
         }
         else if(mid * mid > val) 
         {
               high = mid;
               mid = (low + high)/2;      
          } 
          else 
          {
               low = mid;
               mid = (low + high)/2;
          }
      }
    }
   return 0;
}

$ ./a.out 
Input the values of the number whose square root need be found, and the tolerance limit in the number of digits, of integer variable: tolerance, giving tolerance limit as: 10^{-tolerance}
Enter the number, followed by the # tolerance digits.
2 6
Square root of 2= 1.41421
$ ./a.out 
Input the values of the number whose square root need be found, and the tolerance limit in the number of digits, of integer variable: tolerance, giving tolerance limit as: 10^{-tolerance}
Enter the number, followed by the # tolerance digits.
2 10
Square root of 2= 1.414213562
$ ./a.out 
Input the values of the number whose square root need be found, and the tolerance limit in the number of digits, of integer variable: tolerance, giving tolerance limit as: 10^{-tolerance}
Enter the number, followed by the # tolerance digits.
2 20
^C
$ ./a.out 
Input the values of the number whose square root need be found, and the tolerance limit in the number of digits, of integer variable: tolerance, giving tolerance limit as: 10^{-tolerance}
Enter the number, followed by the # tolerance digits.
2 15
Square root of 2= 1.41421356237309

Beware of entering tolerances finer than the decimal precision of long double.
It may get stuck in a loop that never converges to within the tolerance you request.

> The reason for that is unclear, as know that can get
> infinite number of decimal digits in the decimal expansion of sqrt{2}.
If you're looking for hundreds of digits, then floating point is the wrong tool for the job.

[jonnin]

it may be worth your time to print the size of long double. It varies depending on the underlying hardware and compiler.

[ajiten73]

@salem_c Seems the header file iomanip is the key.
But, your program fails to work (as your version works) even if just add an output file to store output. This is additional to the cout statements.
You can check the working of the below program.

#include <iostream>
#include <iomanip>
#include <fstream>
#include <math.h>
using namespace std;

int main() {
   int val, prec;
   ofstream out_file;
   out_file.open("output_bisect.txt");
   cout << "Enter the number, and precision"<<endl;
   out_file << "Enter the number, and precision"<<endl;
   cin >>val>> prec;
   
   if(val< 0) 
       exit(0);
   else 
   {
      long double low = 0, high = val;
      long double  mid = (low + high)/2;
      int c = 0; long double magn=0;

      while (!c) 
      {
         if((mid*mid-val)<0)
             magn = val - mid*mid;
         
         if (magn <= pow(10, -1*prec))
         {
            cout << "Square root is: " << mid <<endl;
            out_file << "Square root is: " << mid <<endl;
            c = 1;
         }
         else if(mid * mid > val) 
         {
               high = mid;
               mid = (low + high)/2;
          } 
          else 
          {
               low = mid;
               mid = (low + high)/2;  
          }
       }
    }
   return 0;
}

[salem c]

> Seems the header file iomanip is the key.
Yeah, but you still need to call setprecision.

[George P]

The default I/O precision for floating point types is set by the compiler, commonly set to 6 digits. The I/O precision can be be changed by using <iomanip>'s std::setprecision manipulator.

https://en.cppreference.com/w/cpp/io/manip/setprecision

Your OS, the underlying computer hardware and compiler may use a different maximum precision number, the largest machine-native floating point on many platforms is 80 bits.

If you need a greater number of significant digits than your compiler's default there are 3rd party libraries available. GNU's MPFR library, for instance.

https://www.mpfr.org/

Or Boost's multiprecision library. Integer and floating point types available.

https://www.boost.org/doc/libs/1_81_0/libs/multiprecision/doc/html/index.html

[George P]

C++20 has std::format that allows for modifying the output precision for floating point types, the default is set to the maximum allowed number of digits.

#include <iostream>
#include <iomanip>
#include <numbers>
#include <limits>
#include <format>

int main( )
{
   std::cout << "float: " << sizeof(float) << '\n';
   std::cout << std::numeric_limits<float>::digits10 << "\n\n";

   std::cout << "double: " << sizeof(double) << '\n';
   std::cout << std::numeric_limits<double>::digits10 << "\n\n";

   std::cout << "long double: " << sizeof(long double) << '\n';
   std::cout << std::numeric_limits<long double>::digits10 << "\n\n";

   // default precision, outputs 3.14159
   std::cout << std::numbers::pi << '\n';

   // maximum precision, outputs 3.141592653589793
   constexpr auto max_precision { std::numeric_limits<long double>::digits10 + 1 };
   std::cout << std::setprecision(max_precision) << std::numbers::pi << "\n\n";

   // default precision, outputs 3.141592653589793
   std::cout << std::format("{}\n", std::numbers::pi);

   // https://en.cppreference.com/w/cpp/utility/format/formatter#Standard_format_specification
   // outputs 3.141592653589793115997963468544185161590576171875
   // is this really a valid number, though?
   std::cout << std::format("{:2.48f}\n", std::numbers::pi);
}

float: 4
6

double: 8
15

long double: 8
15

3.14159
3.141592653589793

3.141592653589793
3.141592653589793115997963468544185161590576171875

I use MSVC, a long and a long double use the same number of bytes and have the same number of significant digits maximum precision.

If your compiler doesn't support <format> remove lines 25-31 and the <format> header. Compile and run will tell you what your system and compiler can do.

If your compiler doesn't support <format> consider using the {fmt} library, it formed the basis for C++20's stdlib <format> library.

https://fmt.dev/latest/index.html

[ajiten73]

@salim_c: You stated:
>Beware of entering tolerances finer than the decimal precision of long double.
>It may get stuck in a loop that never converges to within the tolerance you >request.

Please elaborate, as have :

if((mid*mid-val)<0)
    magn = val - mid*mid;
         
if (magn <= pow(10, -1*prec))
{

then how having value of: pow(10, -1*prec) less than decimal precision of long double, can cause infinite loop.

[ajiten73]

@George_P Thanks a lot for excellent inputs on the topic, for which even books seems have nothing to impart.
Till now, got inputs from SO (https://stackoverflow.com/q/554063/555045), about an alternative way:
- add in header:
(i) header file : #include <limits>
(ii) typedef std::numeric_limits< double > dbl
- in code, state:
cout.precision(dbl::max_digits10)

This works too, but surprisingly not stated here.
Request if there are gains in one approach over the other; & hence not stated by you?
Also, if outstream's precision is changed, then if store output in file, there the same change should occur too.
But, failed with the SO approach, stated here.

Haven't tried so far your method, but think it would also fail for output stored in file.

Will update soon.

[seeplus]

Although long double has a higher precision than double

This is not mandated by the standard - only that long double has at least the same precision as double. There is no requirement that it has more. For MS VS double and long double both have size of 8 bytes and both the same precision.

[lastchance]

How many digits would you like?

#include <iostream>
#include <vector>
using namespace std;


void fix( vector<int> &V )             // sort out any "carry" operations
{
   for ( int i = 0, carry = 0; i < V.size() || carry > 0; i++ )
   {
      if ( i < V.size() ) carry += V[i];
      int digit = carry % 10;
      carry /= 10;
      if ( i < V.size() ) V[i] = digit;
      else                V.push_back( digit );
   }
}

ostream &operator << ( ostream &out, const vector<int> &V )
{
   out << V.back() << ".";
   for ( int i = V.size() - 2; i >= 0; i-- ) out << V[i];
   return out;
}


int main()
{
   int N = 2;                                               // SINGLE digit number (for the moment)
   int DP = 60;                                             // Decimal places for square root
   int first = 0;
   for ( int digit = 1; digit * digit <= N; digit++ ) first++;  // Find whole-number part of square root
   vector<int> S = { first }, SS = { first * first };       // These hold the prospective square root and its square in REVERSE order

   for ( int decimals = 1; decimals <= DP; decimals++ )     // Add one decimal place at a time
   {
      SS.insert( SS.begin(), { 0, 0 } );                    // The square must add two decimal places
      int nextDigit = 0;
      auto SSnext = SS;
      for ( int digit = 1; digit < 10; digit++ )            // Find the largest final digit such that the square doesn't exceed N
      {
         vector<int> TEST = SS;
         TEST[0] += digit * digit;                                              // These two lines obviate having to multiply 
         for ( int i = 0; i < S.size(); i++ ) TEST[i+1] += 2 * digit * S[i];    //      whole big numbers, which is expensive
         fix( TEST );                                                           // Make sure that all digits are less then 10
         if ( TEST.size() > 2 * decimals + 1 || TEST[2*decimals] >= N ) break;  // Too large
         SSnext = TEST;                                                         // Update the candidate square
         nextDigit = digit;                                                     //      and last digit of the square root
      }
      SS = SSnext;
      S.insert( S.begin(), nextDigit );
      cout << S << '\n';
   }
}

1.4
1.41
1.414
1.4142
1.41421
1.414213
1.4142135
1.41421356
1.414213562
1.4142135623
1.41421356237
1.414213562373
1.4142135623730
1.41421356237309
1.414213562373095
1.4142135623730950
1.41421356237309504
1.414213562373095048
1.4142135623730950488
1.41421356237309504880
1.414213562373095048801
1.4142135623730950488016
1.41421356237309504880168
1.414213562373095048801688
1.4142135623730950488016887
1.41421356237309504880168872
1.414213562373095048801688724
1.4142135623730950488016887242
1.41421356237309504880168872420
1.414213562373095048801688724209
1.4142135623730950488016887242096
1.41421356237309504880168872420969
1.414213562373095048801688724209698
1.4142135623730950488016887242096980
1.41421356237309504880168872420969807
1.414213562373095048801688724209698078
1.4142135623730950488016887242096980785
1.41421356237309504880168872420969807856
1.414213562373095048801688724209698078569
1.4142135623730950488016887242096980785696
1.41421356237309504880168872420969807856967
1.414213562373095048801688724209698078569671
1.4142135623730950488016887242096980785696718
1.41421356237309504880168872420969807856967187
1.414213562373095048801688724209698078569671875
1.4142135623730950488016887242096980785696718753
1.41421356237309504880168872420969807856967187537
1.414213562373095048801688724209698078569671875376
1.4142135623730950488016887242096980785696718753769
1.41421356237309504880168872420969807856967187537694
1.414213562373095048801688724209698078569671875376948
1.4142135623730950488016887242096980785696718753769480
1.41421356237309504880168872420969807856967187537694807
1.414213562373095048801688724209698078569671875376948073
1.4142135623730950488016887242096980785696718753769480731
1.41421356237309504880168872420969807856967187537694807317
1.414213562373095048801688724209698078569671875376948073176
1.4142135623730950488016887242096980785696718753769480731766
1.41421356237309504880168872420969807856967187537694807317667
1.414213562373095048801688724209698078569671875376948073176679

Apparently sqrt(2) to 1000 decimal places is
1.4142135623730950488016887242096980785696718753769480731766797379907324784621070388503875343276415727350138462309122970249248360558507372126441214970999358314132226659275055927557999505011527820605714701095599716059702745345968620147285174186408891986095523292304843087143214508397626036279952514079896872533965463318088296406206152583523950547457502877599617298355752203375318570113543746034084988471603868999706990048150305440277903164542478230684929369186215805784631115966687130130156185689872372352885092648612494977154218334204285686060146824720771435854874155657069677653720226485447015858801620758474922657226002085584466521458398893944370926591800311388246468157082630100594858704003186480342194897278290641045072636881313739855256117322040245091227700226941127573627280495738108967504018369868368450725799364729060762996941380475654823728997180326802474420629269124859052181004459842150591120249441341728531478105803603371077309182869314710171111683916581726889419758716582152128229518488472