RPCSX/llvm
1
0
Fork 0
mirror of https://github.com/RPCSX/llvm.git synced 2024-11-25 20:59:51 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

Use a better algorithm for rounding sqrt results. Change the FIXME about

Browse Source 
this to a NOTE: because pari/gp results start to get rounded incorrectly
after 192 bits of precision. APInt and pari/gp never differ by more than
1, but APInt is more accurate because it does not lose precision after 192
bits as does pari/gp.


git-svn-id: https://llvm.org/svn/llvm-project/llvm/trunk@34834 91177308-0d34-0410-b5e6-96231b3b80d8
This commit is contained in:
Reid Spencer 2007-03-02 04:21:55 +00:00
parent c6a28fcf5d
commit f09aef7698
1 changed files with 12 additions and 8 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

20
lib/Support/APInt.cpp
View File

@ -1239,19 +1239,23 @@ APInt APInt::sqrt() const {
} }
// Make sure we return the closest approximation // Make sure we return the closest approximation
// FIXME: This still has an off-by-one error in it. Test case: // NOTE: The rounding calculation below is correct. It will produce an
// 190 bits: sqrt(694114394047834196220892040454508646882614255319893124270) = // off-by-one discrepancy with results from pari/gp. That discrepancy has been
// 26346050824513229049493703285 (not 26346050824513229049493703284) // determined to be a rounding issue with pari/gp as it begins to use a
// floating point representation after 192 bits. There are no discrepancies
// between this algorithm and pari/gp for bit widths < 192 bits.
APInt square(x_old * x_old); APInt square(x_old * x_old);
APInt nextSquare((x_old + 1) * (x_old +1)); APInt nextSquare((x_old + 1) * (x_old +1));
if (this->ult(square)) if (this->ult(square))
return x_old; return x_old;
else if (this->ule(nextSquare)) else if (this->ule(nextSquare)) {
if ((nextSquare - *this).ult(*this - square)) APInt midpoint((nextSquare - square).udiv(two));
return x_old + 1; APInt offset(*this - square);
else if (offset.ult(midpoint))
return x_old; return x_old;
else else
return x_old + 1;
} else
assert(0 && "Error in APInt::sqrt computation"); assert(0 && "Error in APInt::sqrt computation");
return x_old + 1; return x_old + 1;
} }