numeric_limits<double>::epsilon() gives the error of (double)1.0
in the sense that 1+epsilon is the next double value which is represented.
But how can we calculate the error of arbitrary double number?
numeric_limits<double>::epsilon() gives the error of (double)1.0
in the sense that 1+epsilon is the next double value which is represented.
But how can we calculate the error of arbitrary double number?
Given a number x, you want to find epsilon with the following two properties:
1) x + epsilon != x
2) If 0 <= y < epsilon, x + y == x.
It should not be terribly hard to find the value of epsilon by binary search.
Yes, but its very slow for me.
But an estimate is good for me.
We can use the known representation of double numbers.
If you use the known representation of double numbers, your program will work only on machines that use that representation. If you use the binary-search technique I suggested, it will run anywhere.
Admittedly the binary search is slower. What is it about this particular operation that you need to be able to do quickly?
I want to use error handling, which means every number has a mean value and an absolute error, which can be a real error or error due to the representation.
I think that basically the double representation is the same on every platform. Or the only variable is the number of used bytes. But this is not a serious problem. So, can anyone tell me a function which can be used as an estimate?
I have tried the following code:
#include <iostream>
#include <limits>
using namespace std;
int main()
{
for(double x = 1e-10; x < 1e11; x *= 10)
{
cout << x;
double s = 0;
for(double y = x; y == x; y = x+(s+=x*1e-20)) {}
cout << " " << s/x << endl;
}
cout << "epsilon(): " << numeric_limits<double>::epsilon() << endl;
return 0;
}
The output is:
1e-010 6.46e-017
1e-009 1.034e-016
1e-008 8.28e-017
1e-007 6.63e-017
1e-006 1.059e-016
1e-005 8.48e-017
0.0001 6.78e-017
0.001 1.085e-016
0.01 8.67e-017
0.1 6.95e-017
1 1.111e-016
10 8.89e-017
100 7.11e-017
1000 5.69e-017
10000 9.1e-017
100000 7.28e-017
1e+006 5.83e-017
1e+007 9.32e-017
1e+008 7.46e-017
1e+009 5.97e-017
1e+010 9.55e-017
epsilon(): 2.22045e-016
Why is the calculated value for 1 is only the half of epsilon()?
I think that a good estimation is error x = x*1e-16;