Optimization and double comparison

Optimization and double comparison

[Torquil Macdonald Sørensen]

Hi, is it expected that the attached program only works when no optimization 
is used? I do not understand why the if-test is ever entered in this program. 
The variables that are compared should be exactly the same, since one is used 
to set the value of the other another place in the program. I will attach the 
whole program below.

In short, how can it be that an if-test of the following form can print two 
exactly equal numbers?:

if( a != b) {
	cout << setprecision(70);
	cout << a << " " << b << endl;
} 

Btw, I know that normally one should do a comparison within a certain error 
when comparing doubles, but I am only interested in how the if-test above can 
print two equal numbers?

Here is the program (the if-test is entered when N=4, but not when N=2,3). In 
short, first pot[m][n] is set using calc(m,n), but afterwards the if-test says 
that they are different, if I compile with -O1:

#include <iostream>
#include <iomanip>

using namespace std;

const int M = 2, N = 4;

double calc(const int m, const int n)
{
	double pot = 0.0;
	for(int n2 = 0; n2 < n; ++n2) {
		pot += 0.1;
	}
	return(pot);
}

void check(double pot[][N])
{
	for(int m = 0; m != M; ++m) {
		for(int n = 0; n != N; ++n) {
			cout << "n = " << n << endl;
			if( calc(m,n) != pot[m][n] ) {
				cout << "Error!" << endl;
				cout << setprecision(70);
				cout << pot[m][n] << " and " << calc(m,n) << endl;
			}
		}
	}
}

int main()
{
	double pot[M][N];

	for(int m = 0; m != M; ++m) {
		for(int n = 0; n != N; ++n) {
			pot[m][n] = calc(m,n);
		}
	}

	check(pot);

	return(0);
}

Re: Optimization and double comparison

[Brian Dessent]

Torquil Macdonald Sørensen wrote:

> In short, how can it be that an if-test of the following form can print two
> exactly equal numbers?:
> 
> if( a != b) {
>         cout << setprecision(70);
>         cout << a << " " << b << endl;
> }

You didn't mention what target this is.  If it's x86 then this is the
classic issue of excess precision, wherein operations within the 387
unit occur in extended (80 bit) precision but when transferred out of
the FP unit and stored to memory are truncated to double precision (64
bit.)  If you compare a value that has just been computed with one that
has been previously computed and then stored to memory, they won't
necessarily be equal since one has been truncated while the other still
has its excess precision.  But by passing it on the stack to cout you
essentially force a memory store which discards the excess precision so
they both appear the same.

Much has already been written on this topic, so I suggest just reading
PR323 or googling.  There are numerous workarounds: use sse2 instead of
387, set the 387 to double precision (e.g. -mpc64 or _FP_{GET,SET}CW),
use -ffloat-store, don't compare for absolute equality but against some
small delta, etc.

Brian

Re: Optimization and double comparison

[Torquil Macdonald Sørensen]

On Tuesday 5. August 2008 00:20, Brian Dessent wrote:
> Torquil Macdonald Sørensen wrote:
> > In short, how can it be that an if-test of the following form can print
> > two exactly equal numbers?:
> >
> > if( a != b) {
> >         cout << setprecision(70);
> >         cout << a << " " << b << endl;
> > }
>
> You didn't mention what target this is.  If it's x86 then this is the
> classic issue of excess precision, wherein operations within the 387
> unit occur in extended (80 bit) precision but when transferred out of
> the FP unit and stored to memory are truncated to double precision (64
> bit.)  If you compare a value that has just been computed with one that
> has been previously computed and then stored to memory, they won't
> necessarily be equal since one has been truncated while the other still
> has its excess precision.  But by passing it on the stack to cout you
> essentially force a memory store which discards the excess precision so
> they both appear the same.
>
> Much has already been written on this topic, so I suggest just reading
> PR323 or googling.  There are numerous workarounds: use sse2 instead of
> 387, set the 387 to double precision (e.g. -mpc64 or _FP_{GET,SET}CW),
> use -ffloat-store, don't compare for absolute equality but against some
> small delta, etc.
>
> Brian

Hi Brian, thank you very much. Yes the target is x86. Your answer was very 
interesting. The comparison was just to make sure that a certain quantity 
still had the value that it is supposed to have, so it was only for debugging 
my algorithm. It should be allright if I instead allow a small error in the 
comparison, so I think that I don't have to worry about this problem then.

Thanks again,
Torquil Sørensen

Re: Optimization and double comparison

[Bob Plantz]

On Tue, 2008-08-05 at 10:19 +0200, Torquil Macdonald Sørensen wrote:
> On Tuesday 5. August 2008 00:20, Brian Dessent wrote:
> >
> > Much has already been written on this topic, so I suggest just reading
> > PR323 or googling.  There are numerous workarounds: use sse2 instead of
> > 387, set the 387 to double precision (e.g. -mpc64 or _FP_{GET,SET}CW),
> > use -ffloat-store, don't compare for absolute equality but against some
> > small delta, etc.
> >
> > Brian
> 
> Hi Brian, thank you very much. Yes the target is x86. Your answer was very 
> interesting. The comparison was just to make sure that a certain quantity 
> still had the value that it is supposed to have, so it was only for debugging 
> my algorithm. It should be allright if I instead allow a small error in the 
> comparison, so I think that I don't have to worry about this problem then.

You may wish to at least note that the default for 32-bit compile is
x87, but for 64-bit it is SSE2. As Brian noted, x87 uses 80-bit
internally, but SSE2 uses only 64-bit.

Somebody down the road may port this to 64-bit and get different
results. (In my experience, the "somebody down the road" may be myself,
and I've forgotten why I did something. :) )

Bob