C++Talk.NET

archondas · Guest

I have encountered a very strange behaviour
of floating point accuracy depending on the way
I allocate memory! The results are slightly different
but should they? I used the following simple program:

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

#define STATIC_SIZE 60000000
double static_array[STATIC_SIZE];

int main(int argc, char* argv[])
{
if (argc < 2)
{
cout << "usage: " << argv[0] << " nn";
cout << "example: " << argv[0] << " 10000000n";
exit(1);
}

int n = (int) strtol(argv[1], NULL, 0);

double sum = 0.0;
double static_sum = 0.0;
double dynamic_sum = 0.0;

double* dynamic_array = new double[n];
const double pi = 4.0*atan(double(1.0));

for (int i = 0; i < n; i++)
{
double value_i = cos(pi*i/(n-1));
dynamic_array[i] = value_i;
static_array[i] = value_i;
dynamic_sum += dynamic_array[i];
static_sum += static_array[i];
sum += value_i;
}

cout.setf(ios::scientific, ios::floatfield);
cout.precision(16);
cout << "sum = " << sum << endl;
cout << "static_sum = " << static_sum << endl;
cout << "dynamic_sum = " << dynamic_sum << endl;

delete[] dynamic_array;
return 0;
}

Here is the output with g++-3.4.3 -O2
--------------------------------------
# ./run 100000
sum = 1.0667244865203429e-11
static_sum = 1.0667244865203429e-11
dynamic_sum = -3.8704595084482207e-12

# ./run 1000000
sum = -4.6566472899911560e-10
static_sum = -4.6566472899911560e-10
dynamic_sum = -3.4879821253497312e-10

# ./run 10000000
sum = 8.2471927065341788e-09
static_sum = 8.2471927065341788e-09
dynamic_sum = 8.3632101244290880e-09

# ./run 40000000
sum = -1.0021416363414914e-07
static_sum = -1.0021416363414914e-07
dynamic_sum = -1.0021962060235978e-07

exactly the same results are produced with
icc -O2 (icc version 8.1.026).

For the test I used an AMD ATHLON-XP.

Any comments?

[ See http://www.gotw.ca/resources/clcm.htm for info about ]
[ comp.lang.c++.moderated. First time posters: Do this! ]

Naveen · Guest

Hi archondas.
I think the discrepancy might be due to some compiler optimization. I
see that you have compiled with O2 flag.

Try compiling without the O2 flag. The results should be the same.
--
Regards
Naveen.

[ See http://www.gotw.ca/resources/clcm.htm for info about ]
[ comp.lang.c++.moderated. First time posters: Do this! ]

Ron Natalie · Guest

archondas wrote:

Ulrich Eckhardt · Guest

archondas wrote:

archondas · Guest

Ok for all of you,
here is the new code! I think it is more fair than
the previous one! Take a look!

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

#define STATIC_SIZE 60000000
double static_array[STATIC_SIZE];

int main(int argc, char* argv[])
{
if (argc < 2)
{
cout << "usage: " << argv[0] << " nn";
cout << "example: " << argv[0] << " 10000000n";
exit(1);
}

int n = (int) strtol(argv[1], NULL, 0);

double sum = 0.0;
double static_sum = 0.0;
double dynamic_sum = 0.0;

double* dynamic_array = new double[n];
const double pi = 4.0*atan(double(1.0));

for (int i = 0; i < n; i++)
{
dynamic_array[i] = cos(pi*i/(n-1));
static_array[i] = cos(pi*i/(n-1));
}

for (int i = 0; i < n; i++)
{
double value_i = cos(pi*i/(n-1));
dynamic_sum += dynamic_array[i];
static_sum += static_array[i];
sum += value_i;
}

cout.setf(ios::scientific, ios::floatfield);
cout.precision(16);
cout << "sum = " << sum << endl;
cout << "static_sum = " << static_sum << endl;
cout << "dynamic_sum = " << dynamic_sum << endl;

delete[] dynamic_array;
return 0;
}

# ./run 1000
sum = 8.4376949871511897e-14
static_sum = 8.3488771451811772e-14
dynamic_sum = 8.3488771451811772e-14

# ./run 10000
sum = -9.6567198681896116e-13
static_sum = -5.5733195836182858e-14
dynamic_sum = -5.5733195836182858e-14

# ./run 100000
sum = 1.0667244865203429e-11
static_sum = -3.8704595084482207e-12
dynamic_sum = -3.8704595084482207e-12

# ./run 1000000
sum = -4.6566472899911560e-10
static_sum = -3.4879821253497312e-10
dynamic_sum = -3.4879821253497312e-10

# ./run 10000000
sum = 8.2471927065341788e-09
static_sum = 8.3632101244290880e-09
dynamic_sum = 8.3632101244290880e-09

# ./run 40000000
sum = -1.0021416363414914e-07
static_sum = -1.0021962060235978e-07
dynamic_sum = -1.0021962060235978e-07

Without -O2 for g++ the results are the same for all the
above cases. But as you see for very large numbers
the result is closer to zero for sum than for static_sum
or dynamic_sum. Should it? Yes of course this is
closer to the correct answer which is zero! But why isn't it
for smaller numbers? Why floating point operations loose
accuracy when storing the result in a local variable for
indices (n <10000000) and not for (n>=10000000). What changes? And guys,
why should the result depend on -O2.
Why optimizing the code should produce different results!
Let's try long double instead of double.
That worked fine! Both with or without -O2.
What is the difference guys?
icc with long double workder fine even with -O2.
icc with double behaves just like gcc giving the same
results.

Summary:
--------
Floating point addition of #n doubles is more accurate
when the doubles are read from memory for n<10000000.
Less accurate when n >= 10000000. (with -O2 enabled)

With -O2 disabled there is not difference although
the results seem to be less accurate than with -O2 enabled.

long doubles do not have such problems.

Conclusion:
-----------
Some guys should fix their -O2 implementation.
Use long doubles to be on the safe side!

[ See http://www.gotw.ca/resources/clcm.htm for info about ]
[ comp.lang.c++.moderated. First time posters: Do this! ]

Jeroen · Guest

This definitely has to do with the optimizations (-O2 switch) you are
using. If you compile without optimizations, you'll see there is no
difference between the three results at all, but also that the result
is different than ALL three results you get when optimizations are
turned on.

Floating point numbers have limited precision. Combined with
optimizations, this can lead to small differences in results. What
you're seeing here is a very tiny difference in the value of value_i
(Order 10^-17), accumulated over thousands of iterations.

[ See http://www.gotw.ca/resources/clcm.htm for info about ]
[ comp.lang.c++.moderated. First time posters: Do this! ]

Antoun Kanawati · Guest

archondas wrote:

pierrebai@hotmail.com · Guest

It is not possible to be certain without looking at the disassembly,
but I would suspect that the static array value is not reloaded after
the store while the dynamic value is. Must be optimized because the
static array is visible.

[ See http://www.gotw.ca/resources/clcm.htm for info about ]
[ comp.lang.c++.moderated. First time posters: Do this! ]

Ron Natalie · Guest

archondas wrote:

archondas · Guest

Dear Ron,
I do not know from which scientific area you come
from, but have in mind that the slightest change
in floating point numbers can have a great impact
when somebody is trying to solve a very ill-conditioned
linear system. Such systems arise from direct solution
methods applied in scattering problems.

So as for the million floating point values, let
me ask you something. Should it make a difference
when adding from a static allocated vector or from
a dynamic allocated one? Should those result be different
when adding directly to a variable? I think not!
Adding millions of values should only make a difference
for the outcome of the addition. In our case the result
is zero but the outcome of the addition is not zero.
This type of error has to do with the 64-80 bit of FPU
registers. But when you are using a specific type
double, long double, then, each floating point
truncation error when converting cos(something) to
a number is always the same and depends on the value of
cos(something). Having that in mind the results from
our example should always be the same! Why they are not
when we are using type double and -O2 switch is something
we are investigating here! As you remember -O switch gives
as everywhere the same results. What changes with -O2?

You said
--------

Ulrich Eckhardt · Guest

archondas wrote:

Ron Natalie · Guest

archondas wrote:

archondas · Guest

Dear Ron,
Now it is clear why you cannot understand the importance
of the accuracy in floating point calculations! Engineers
usually are satisfied if they have achieved the solution
with an error of 10%. In my case in order to achieve
acceptable results then you have to solve the linear system
as much accurately as possible. So I think that it is clear why changes in
the results like those I described above are not acceptable!

On the other hand here we are examining whether memory
allocation in one way or another should influence
or not the results and under which circumstances the
results are as much accurate as possible(register variables, memory
variables ...). So discussing whether millimiters
make a difference when we are measuring kilometers is neither C++ nor a
compiler issue! You should address this
in some other (list or forum).

[ See http://www.gotw.ca/resources/clcm.htm for info about ]
[ comp.lang.c++.moderated. First time posters: Do this! ]

Francis Glassborow · Guest

In article
<17606844081197021a9833e0822d3c3e (AT) localhost (DOT) talkaboutprogramming.com>,
archondas <dkouroun (AT) cc (DOT) uoi.gr> writes

kanze@gabi-soft.fr · Guest

archondas wrote: