Groups | Search | Server Info | Keyboard shortcuts | Login | Register [http] [https] [nntp] [nntps]

Groups > comp.lang.c > #138984

Re: storing and retrieving a double may show a difference?

Show all headers | View raw

On 12/22/2018 10:08 PM, Ben Bacarisse wrote:
> Manfred <noname@add.invalid> writes:
> 
>> On 12/21/2018 6:10 PM, Ben Bacarisse wrote:
<snip>
>>> No, it's a straight forward logic error.  You change A[1][1] after
>>> calculating diagonal[1] so when you compare the two, they are no longer
>>> the same.
>> Besides, also by simple code inspection, there appears to be more to
>> be corrected - error propagation is not correctly estimated, which
>> affects the assertions. More precisely at least two ingredients of
>> this are missing:
> 
> Yes, but I think the asserts are there only because an error was
> spotted.  I'd go further -- it's rarely useful to use any floating point
> comparisons in an assert because picking the boundary between OK and not
> OK is almost impossible.
I don't fully agree. I find asserts in floating point math very useful, 
but they must be used in the proper way, which is (the proper way) 
somewhat counter-intuitive in that you must assume there are no /exact/ 
values, which is odd in a binary world.
It is true that picking the exact boundary is quite hard, but still some 
reasonable estimate of rounding errors is possible, and with such, 
asserts can help spot errors - you can't get 100% error detection (but 
no test does), but you can get a significant rate.

   Much better to do a full error analysis or to
> collect error data as you go and present that as part of the output.
This can be an additional debug strategy, yes.

   Of
> course, if this is some sort of real-time calculation for a game, say,
> then you would probably not do either.
> 
>> - using DBL_EPSILON
>> - considering how many operations are involved in the checked result
>> that can affect accuracy (esp. additions and subtractions)
>>

<snip>
>> /****
>> Ben is right that there is a logic error, but still in floating point == and != tests require extra attention: in fact they are more correctly used as /identity/ tests rather than /compare equal/ tests.
>> Here this is correct only for architectures that are guaranteed to
>> preserve the FP binary representation.
> 
> I'm not entirely sure I know what you mean here.  In C, == and != /are/
> equality tests and they should ignore any changes or differences in
> representation.  For example, negative zero should equal positive zero,
> and 1, say, in any bit-length should equal 1 in any other bit-length.
You are technically right about this, but my point is about practical 
use of == and != .
In practice, since FP numbers are typically (at least all that need be 
verified) the result of some expression, and they are affected by 
rounding errors, it makes no sense to use == to compare the result of 
two expressions that /in theory/ do compare equal, but in finite 
precision we know they don't.
We know we must use something like (fabs(x-v) < something*DBL_EPSILON), 
and in this perspective the equivalence of positive zero and negative 
zero is not really relevant.
My point about using == as an identity test above is substantially about 
using it to compare a number with itself.
In the example at hand, this means verifying e.g. that an index in an 
array of FP numbers points to the expected element, or that some value 
is not modified between the point it is stored and the point it is 
retrieved (which is the problem you spotted).
The identity test, as a further remark, is guaranteed to work only on 
architectures that do not have the 80/64 bit issue, or something else 
that does not preserve the binary representation of a number between 
storing and retrieving.

> 
> I think I know what you are getting at, but I don't think you've really
> explained it, and I won't even try since these are complex matters and I
> don't have anything like enough time!
Yes, I did not explain it, I only recalled the conclusion. Done now, in 
case anyone is interested.

> 
>> ****/
>>>> 		if (i > 0 && x != diagonal[i]){
>>>> 			printf("i %d, x %.15g != diagonal %.15g\n",i, x ,diagonal[i]);
>>>> 			assert(0);
>>>> 		}
>>>> 		for (j = i; j < height_A; ++j) {
>>>> 			*(A + j * width_A + i) /= x;
>>>> 		}
>>>> 		assert(*(A + i * width_A + i) == 1.0);
>>>> 		// now diag elem = 1;
>>>> 		for (k = i + 1; k < width_A; ++k) {
>>>> 			x = *(A + i * width_A + k); // mult
>>>> 			for (j = i; j < height_A; ++j) {
>>>> 				*(A + j * width_A + k) -= x * *(A + j * width_A + i);
>>>
>>> **** Here you can have k == 1 and j == 1 and you change the element of A
>>> that is on the diagonal.
>>>
>>>> 			}
>> /*****************
>> This is another example:
>> - The constant 1.0e-300 is arbitrary. The test should be based on DBL_EPSILON an (some estimate of) the size of x.
>> - Is this meant to be a test like == 0 or < 0 ? If it is the former, then a fabs() call is missing
>> ******************/
>>>> 			assert(*(A + i * width_A + k) < 1.0e-300);
> 
> I agree about the test, but I'd say get rid of the assert!  Either ignore
> small non-zero values, do a more thorough analysis, or accumulate data
> on the errors for later reporting.
> 
> Only the middle option is really robust.  There are a ill-conditioned
> matrices that produce entirely unreliable results using this method, but
> there are also well-understood tests for these.

Maybe we see a different purpose for the assert. I think it is 
legitimate to use it in places like this, where it could spot some 
logical errors.
OTOH testing for /accuracy/ of the result, like "is this algorithm 
suitable for the kind of matrices it is used with?" is a different 
matter, and you are right that a simple assert() can't do, while the 
alternative you indicate is more appropriate.

> 
>>>> 		}
>>>> 	}
>>>
>>> <snip rest of function>
>>>> 	return x;
>>>> }
>>>
>>
> 

Back to comp.lang.c | Previous | Next — Previous in thread | Next in thread | Find similar | Unroll thread

Thread

storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-20 23:25 -0800
  Re: storing and retrieving a double may show a difference? fir <profesor.fir@gmail.com> - 2018-12-21 00:56 -0800
    Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-21 07:02 -0800
      Re: storing and retrieving a double may show a difference? Bart <bc@freeuk.com> - 2018-12-21 15:33 +0000
        Re: storing and retrieving a double may show a difference? "Rick C. Hodgin" <rick.c.hodgin@gmail.com> - 2018-12-21 11:16 -0500
          Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-21 09:05 -0800
            Re: storing and retrieving a double may show a difference? fir <profesor.fir@gmail.com> - 2018-12-21 09:43 -0800
            Re: storing and retrieving a double may show a difference? "Rick C. Hodgin" <rick.c.hodgin@gmail.com> - 2018-12-21 12:57 -0500
            Re: storing and retrieving a double may show a difference? Ben Bacarisse <ben.usenet@bsb.me.uk> - 2018-12-21 20:06 +0000
      Re: storing and retrieving a double may show a difference? Ben Bacarisse <ben.usenet@bsb.me.uk> - 2018-12-21 17:10 +0000
        Re: storing and retrieving a double may show a difference? fir <profesor.fir@gmail.com> - 2018-12-21 09:17 -0800
          Re: storing and retrieving a double may show a difference? fir <profesor.fir@gmail.com> - 2018-12-21 09:20 -0800
          Re: storing and retrieving a double may show a difference? fir <profesor.fir@gmail.com> - 2018-12-21 09:24 -0800
        Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-21 21:46 -0800
        Re: storing and retrieving a double may show a difference? Manfred <noname@add.invalid> - 2018-12-22 21:41 +0100
          Re: storing and retrieving a double may show a difference? Ben Bacarisse <ben.usenet@bsb.me.uk> - 2018-12-22 21:08 +0000
            Re: storing and retrieving a double may show a difference? Manfred <noname@add.invalid> - 2018-12-23 20:34 +0100
              Re: storing and retrieving a double may show a difference? Tim Rentsch <txr@alumni.caltech.edu> - 2018-12-28 09:14 -0800
          Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-22 19:16 -0800
      Re: storing and retrieving a double may show a difference? fir <profesor.fir@gmail.com> - 2018-12-21 09:12 -0800
      Re: storing and retrieving a double may show a difference? James Kuyper <jameskuyper@alumni.caltech.edu> - 2018-12-21 09:23 -0800
        Re: storing and retrieving a double may show a difference? Ben Bacarisse <ben.usenet@bsb.me.uk> - 2018-12-21 20:03 +0000
          Re: storing and retrieving a double may show a difference? James Kuyper <jameskuyper@alumni.caltech.edu> - 2018-12-22 00:10 -0500
            Re: storing and retrieving a double may show a difference? Ben Bacarisse <ben.usenet@bsb.me.uk> - 2018-12-22 11:06 +0000
  Re: storing and retrieving a double may show a difference? mark.bluemel@gmail.com - 2018-12-21 01:00 -0800
    What type of problem is it? (Was: storing and retrieving a double may show a difference?) gazelle@shell.xmission.com (Kenny McCormack) - 2018-12-21 09:10 +0000
      Re: What type of problem is it? (Was: storing and retrieving a double may show a difference?) mark.bluemel@gmail.com - 2018-12-21 01:33 -0800
  Re: storing and retrieving a double may show a difference? Bart <bc@freeuk.com> - 2018-12-21 10:16 +0000
  Re: storing and retrieving a double may show a difference? Ben Bacarisse <ben.usenet@bsb.me.uk> - 2018-12-21 11:25 +0000
  Re: storing and retrieving a double may show a difference? James Kuyper <jameskuyper@alumni.caltech.edu> - 2018-12-21 08:06 -0500
    Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-21 07:13 -0800
      Re: storing and retrieving a double may show a difference? James Kuyper <jameskuyper@alumni.caltech.edu> - 2018-12-21 09:10 -0800
    Re: storing and retrieving a double may show a difference? Tim Rentsch <txr@alumni.caltech.edu> - 2018-12-28 08:16 -0800
  Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-21 08:09 -0800
  Re: storing and retrieving a double may show a difference? "Chris M. Thomasson" <invalid_chris_thomasson@invalid.invalid> - 2018-12-22 17:23 -0800
    Re: storing and retrieving a double may show a difference? "Chris M. Thomasson" <invalid_chris_thomasson@invalid.invalid> - 2018-12-22 17:58 -0800
  Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-22 19:30 -0800
    Re: storing and retrieving a double may show a difference? Richard Damon <Richard@Damon-Family.org> - 2018-12-22 22:48 -0500
      Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-22 20:10 -0800
        Re: storing and retrieving a double may show a difference? "Chris M. Thomasson" <invalid_chris_thomasson@invalid.invalid> - 2018-12-22 20:41 -0800
        Re: storing and retrieving a double may show a difference? Barry Schwarz <schwarzb@dqel.com> - 2018-12-22 22:02 -0800
          Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-23 02:43 -0800
            Re: storing and retrieving a double may show a difference? Ben Bacarisse <ben.usenet@bsb.me.uk> - 2018-12-23 12:38 +0000
              Re: storing and retrieving a double may show a difference? chanrasjid@gmail.com - 2018-12-23 05:12 -0800
                Re: storing and retrieving a double may show a difference? fir <profesor.fir@gmail.com> - 2018-12-23 05:32 -0800
                Re: storing and retrieving a double may show a difference? Ben Bacarisse <ben.usenet@bsb.me.uk> - 2018-12-23 14:55 +0000
            Re: storing and retrieving a double may show a difference? Richard Damon <Richard@Damon-Family.org> - 2018-12-23 07:43 -0500
    Re: storing and retrieving a double may show a difference? Manfred <noname@add.invalid> - 2018-12-23 20:42 +0100

csiph-web