Thread: Sqrt() innacurate?

  1. #1
    Bored Programmer
    Join Date
    Jul 2009
    Location
    Tomball, TX
    Posts
    428

    Sqrt() innacurate?

    Hey guys,
    I was trying to help someone with a problem involving finding if a number is a perfect square. I came to the conclusion that
    Code:
    int root;
    int square; 
    for(square = 0; square < 100; square++)
    {
      root = sqrt(square);
      if(root * root == square)
      {
        cout<<square;
      }
    }
    Given that square can have any beginning point and ending point within the numeric_limits. His teacher then responded. Due to truncation of sqrt() when dealing with alot of digits you will "miss" some perfect squares. This didn't make sense to me. Any number that requires a truncation of sqrt() wouldn't be a perfect square. To try and demonstrate the error that the teacher said would occur I coded a small program.

    Code:
    #include <iostream>
    #include <stdio.h>
    #include <limits>
    #include <math.h>
    using namespace std;
    
    double y;
    int rootfound;
    int main(void)
    {
      // start at the max value of an integer and find the first square using my method
      for(int x = numeric_limits<int>::max(); x > 0; x--)
      {
        int z = (int)sqrt(x);
        if(z*z == x)
        {
          rootfound = z;
          cout<<"found a perfect square: "<<x<<endl;
          cout<<"its root is "<<z<<endl;
          x = 0;
        } 
      }
      // now test the next perfect square and see if it is over the numeric limits of
      // an integar
      rootfound++;
      cout<<"testing the next number up from this root: "<<endl;
      cout<<"The new root to work with is "<<rootfound<<endl;
      cout<<"Its perfect square is "<<(rootfound*rootfound)<<endl;
      if((rootfound*rootfound) < 0) 
      {
        cout<<"this square was over the numeric limits so we found the first square
                     root lower than limits"<<endl;
      }
      else
        cout<<"this square was under the number limits so the sqrt() function
                    caused us to miss a perfect square"<<endl;
     
      system("pause");
    }
    As I expected the result still finds the first perfect square counting down from the numberic_limits of a int. Can anyone demonstrate how the method used would exclude a perfect square that is within the limits of an int? Further why would sqrt() be in a standard library math.h if it wasn't dependable?
    Last edited by Lesshardtofind; 09-27-2010 at 08:25 AM.
    Virtual reality hello world http://www.rodneybrothers.com/vr/vrh...rld/index.html in html and javascript.
    Viewable with dodocase, google cardboard, OR, and other compatible VR gear.

  2. #2
    The larch
    Join Date
    May 2006
    Posts
    3,573
    I guess due to floating point inaccuracies, it might not be unheard of for the result being slightly lower, e.g sqrt(9) coming out as 2.999999981.

    However, this test program shows that it doesn't happen (at least for me) for the int range.

    Code:
    #include <iostream>
    #include <limits>
    #include <math.h>
    using namespace std;
    
    typedef int IntType ;
    int main(void)
    {
        for (IntType i = 0; i <= (IntType)sqrt(std::numeric_limits<IntType>::max()); ++i) {
            IntType result = (IntType)sqrt(i * i);
            if (result != i) std::cout << i << ' ' << result << '\n';
        }
        std::cout << "All\n";
    }
    I might be wrong.

    Thank you, anon. You sure know how to recognize different types of trees from quite a long way away.
    Quoted more than 1000 times (I hope).

  3. #3
    Bored Programmer
    Join Date
    Jul 2009
    Location
    Tomball, TX
    Posts
    428
    Exactly my point. I can't demonstrate it happening within the integar range on a number that is a perfect square. If it can't be demonstrated I dont' understand a teacher imposing the restrictions, but obviously if its his course you have to follow his rules to pass. It just made me interested in trying to find the error to which he is refering. Maybe its not an accurate error rate like only comes back wrong 1/ 1000 times or someting, but in that case doesn't every programming method have a "possible" error rate due to data corruption and the 1/1,000,000 chance that your computer will misfire.

    Further any other check for a perfect square not using a sqrt() would seem to become time consuming on the CPU running extra loops ect. What if you were dealing with root numbers 0 - 4000. The only other method for approaching I can think of would be
    Code:
    for(square = 0; square < 4000; square++)
    {
      if(IntToCompare == square * square)
      {
        cout<<IntToCompare<<" is a perfect square.";
      }
    }
    Just checking this would be cumbersome. Then add the possibility that you are checking if a list of 2000 integars are a perfect square it seems this would take WAY to long in comparison to the method mentioned before.
    Last edited by Lesshardtofind; 09-27-2010 at 08:24 AM.
    Virtual reality hello world http://www.rodneybrothers.com/vr/vrh...rld/index.html in html and javascript.
    Viewable with dodocase, google cardboard, OR, and other compatible VR gear.

  4. #4
    The larch
    Join Date
    May 2006
    Posts
    3,573
    You could also write your own sqrt (approximation - int precision is enough) which is guaranteed to suit the algorithm.

    Or perform a binary search on the range.
    I might be wrong.

    Thank you, anon. You sure know how to recognize different types of trees from quite a long way away.
    Quoted more than 1000 times (I hope).

  5. #5
    Algorithm Dissector iMalc's Avatar
    Join Date
    Dec 2005
    Location
    New Zealand
    Posts
    6,318
    Actually, if you use floats then you can run into a problem even without using sqrt.
    The following fails when i is 4097 because 16785409 can't be represented as a float.
    It instead stores 16785408.f into f:
    Code:
    	for (unsigned int i=0; i<65536; ++i)
    	{
    		float f = (float)i*i;
    		ASSERT((unsigned int)f == i*i);
    	}
    That is possibly the kind of thing the teacher was worried about.

    Neverthelsss, he/she is perfectly correct that using sqrt is not the ideal way to do this, as that usually entails use of floating point data types, and the issues they bring (performance, accuracy, etc).
    You've recognised the need for a cast of the return value of sqrt, but you also need to recognise that the warning it otherwise produces is not merely there to get in your way, and that it really is trying to tell you something that's often worth paying attention to.
    Last edited by iMalc; 09-28-2010 at 01:48 AM.
    My homepage
    Advice: Take only as directed - If symptoms persist, please see your debugger

    Linus Torvalds: "But it clearly is the only right way. The fact that everybody else does it some other way only means that they are wrong"

  6. #6
    Officially An Architect brewbuck's Avatar
    Join Date
    Mar 2007
    Location
    Portland, OR
    Posts
    7,396
    What makes you think that a floating point can represent any arbitrary integer precisely?

    A double is 64 bits but its range is much, much larger than the range of a 64-bit integer. Therefore, there must be some integers it can't represent.

    Right?
    Code:
    //try
    //{
    	if (a) do { f( b); } while(1);
    	else   do { f(!b); } while(1);
    //}

  7. #7
    Guest Sebastiani's Avatar
    Join Date
    Aug 2001
    Location
    Waterloo, Texas
    Posts
    5,708
    Quote Originally Posted by anon View Post
    You could also write your own sqrt (approximation - int precision is enough) which is guaranteed to suit the algorithm.

    Or perform a binary search on the range.
    Yep, I've found that method to be effective as well. Another thing you can do is optimize the search range, such that it is already quite close to the actual square root. This implementation calculates the largest power of two needed to contain the root, eg: if N = 1234 then H (high) = 128... since the actual root is 111, H is obviously a pretty close initial guess. It might be improved upon by working out the minimum lower-bound as well, although I'm not sure that it's really worth the effort! It seems to run faster than the Newton-Raphson, AFAICT.

    Code:
    size_t isqrt( size_t val )
    {
    	size_t const
    		one = 1;
    	size_t
    		div,
    		res,
    		low = one,
    		hig = one, 
    		tmp = val;	
    	if( !val )
    		return val;
    	while( tmp >>= one )
    		++hig;
    	hig = one << ( ( hig + ( hig & one ) ) >> one );
    	for( ;; )
    	{
    		res = low + ( ( hig - low ) >> one );
    		div = val / res;
    		if( div < res )
    		{
    			if( hig == res )
    				break;
    			hig = res;
    		}
    		else if( div > res )
    		{
    			if( low == res )
    				break;
    			low = res;
    		}
    		else
    			break;
    	}
    	return res;
    }
    Last edited by Sebastiani; 09-28-2010 at 01:32 PM. Reason: simplification
    Code:
    #include <cmath>
    #include <complex>
    bool euler_flip(bool value)
    {
        return std::pow
        (
            std::complex<float>(std::exp(1.0)), 
            std::complex<float>(0, 1) 
            * std::complex<float>(std::atan(1.0)
            *(1 << (value + 2)))
        ).real() < 0;
    }

Popular pages Recent additions subscribe to a feed

Similar Threads

  1. undefined reference to sqrt
    By shamma in forum C Programming
    Replies: 1
    Last Post: 05-09-2009, 01:42 PM
  2. Replies: 5
    Last Post: 06-01-2006, 04:37 PM
  3. sqrt() function help
    By willc0de4food in forum C Programming
    Replies: 5
    Last Post: 03-14-2005, 09:07 PM
  4. SQRT Mystery...
    By KneeLess in forum C Programming
    Replies: 7
    Last Post: 03-23-2004, 07:49 AM
  5. how sqrt ?
    By sambs1978 in forum C Programming
    Replies: 3
    Last Post: 09-20-2001, 08:14 AM