Some dodgy #ifdef and #pragma logic at check-in 5d9e9364808793d6 seems to work around the problem. Please speak up if you see a better solution.

GCC devs report that this issue is a dup of a prior bug. They recommend using either -fexcess-precision=standard or -std=c99 to work around. Those both work for me when entered on the command-line. However, I could not get them to work with #pragma. And we do not want "sqlite3.c" to require special compiler options - it should "just work". Hence, I'll stay with the dodgy #ifdef and #pragma solution unless and until somebody can suggest a better solution.

The history on that "bug", spanning over 23 years now, shows that it is not considered a bug. FWIW, I agree with the analysis and conclusion.

An interesting work-around appears in comment #87, where the use of the "volatile" storage class is used to achieve the same effect as that #pragma, except in a more targeted way.

As I understand the issue, there is little reason to expect that it will be unique to gcc. For that reason, I think judicious use of volatile would be a better fix because it would apply to present and future compilers other than gcc whose optimizers might induce the same behavior.

That bug is 23 years old.... Impressive stuff.

You try using "volatile". I couldn't get it to do anything. Let me know if you find a solution using volatile.

Using your post #1 code as a baseline, I was able to get the "suspect" and "correct" results to match thusly (changes per diff):

3,4c3,4
< void dekkerMul2(double *x, double y, double yy){
<   double hx, tx, hy, ty, p, q, c, cc;
---
> void dekkerMul2(volatile double *x, double y, double yy){
>   volatile double hx, tx, hy, ty, p, q, c, cc;

I found that both the argument and local "volatile" qualifier were needed, and did not try anything more restricted.

This was using gcc 11.3 .

Very good. I didn't try adding volatile to the function parameters. I reproduce your results, except that I omitted volatile from hx and hy.

On the other hand, there is no performance difference. But the executable size is around a hundred bytes larger. So you can argue that "volatile" is "cleaner" and less likely to cause problems on other compilers. But I note:

• The problem manifests as a very small error in binary-to-decimal conversion - literally an error in the 16-th significant digit.

• The problem only comes up on i686. That's pretty much a legacy platform at this point, isn't it?

So I'm not sure yet that I want to make the change. Convince me.

Apparently when compiling for 64-bit, GCC by default emits SSE instrucions but for 32-bit it by default emits 387 instructions. 387 uses extended precision by default (it can be changed by setting the Precision Control Field in the x86 FPU Control Word). If -ffloat-store is enabled then results are always stored into memory and load back for subsequent calculations instead of just using value in the register. SSE uses single or double precision depending on machine instruction opcode. You can choose whether you want GCC to emit SSE or 387 instructions, either by adding -mfpmath=sse or -mfpmath=387 respectively to the command line, or by using pragma:

#pragma GCC push_options
#pragma GCC target("fpmath=sse")

 /* your functions here */

#pragma GCC pop_options

Convince me.

I thought the point regarding compiler independence, together with the (as-yet undemonstrated) notion that the same issue could arise with other compilers, was pretty convincing.

The fact is that the code is doing things where the extra few bits of precision kept in FP registers instead of double memory slots makes a difference. The "volatile" qualifier can be thought of as a way to tell any compiler "The significant part of this quantity is just what is stored in such-and-such memory slot; please use only that in further calculation."

It is much like "#pragma GCC optimize("float-store")", except written in C rather than gcc-speak.

Sorry if I'm just repeating myself.

The "volatile" qualifier can be thought of as a way to tell any compiler "The significant part of this quantity is just what is stored in such-and-such memory slot; please use only that in further calculation

I'm no expert, but that's NOT how I'd try to explain volatile...

The way I think about it instead, is to inform the compiler that a particular variable
(i.e. memory slot) may be modified in ways outside the program's control
(meaning by other means than the code the compiler see, and thus can reason about),
and thus the compiler should refrain from aggresive optimizations on that variable,
which typically means avoiding registers and prefer reloading it from memory instead.

But in practice, I know I don't know enough, and won't rely on volatile :)
I think it's useful/essential to multi-threading experts. Not regular devs. My $0.02.

NEVER use volatile for multi-threading purpose, it does not provide any synchronisation, nor does it ensure the order of execution of operations. volatile was specifically intended to be used when interfacing with memory-mapped hardware, signal handlers, and the setjmp machine code instruction.

FWIW, I didn't mean volatile replaced atomics, but I've seen I believe usage of it in that MT context. In any case, OT, and since SQLite is C89-compatible and the memory model of C wasn't specific until C11, maybe not even portable.

[Grammo fixed via edit.]

I was not attempting an explanation of volatile generally. It might have been clearer for me to write, "The use of the 'volatile' qualifier here can be thought of as ...", or maybe not. The context makes this evident, and needless words can detract from the point being made.

The effect that I claim is one which necessarily flows from the semantics of volatile. After the value is written to the (volatile) memory slot, the code generator is obliged to read it from that slot when the value is later reused via the slot's name. I saw no need to tutor Richard upon this fact. So my sub-point was just a reminder.

Not every architecture/compiler/compile options need volatile here to give correct results, it would be nice to be able to disable use of volatile.

#ifndef SQLITE_FP_USE_VOLATILE
/* Use of voltile enabled by default,
user can define SQLITE_FP_USE_VOLATILE as 0 to disable it. */
#define SQLITE_FP_USE_VOLATILE 1
#endif

#if SQLITE_FP_USE_VOLATILE
#define FP_VOLATILE volatile
#else
#define FP_VOLATILE
#endif

Index: src/func.c
==================================================================
--- src/func.c
+++ src/func.c
@@ -1685,15 +1685,15 @@
 **
 ** Variables are marked "volatile" to defeat c89 x86 floating point
 ** optimizations can mess up this algorithm.
 */
 static void kahanBabuskaNeumaierStep(
-  volatile SumCtx *pSum,
-  volatile double r
+  FP_VOLATILE SumCtx *pSum,
+  FP_VOLATILE double r
 ){
-  volatile double s = pSum->rSum;
-  volatile double t = s + r;
+  FP_VOLATILE double s = pSum->rSum;
+  FP_VOLATILE double t = s + r;
   if( fabs(s) > fabs(r) ){
     pSum->rErr += (s - t) + r;
   }else{
     pSum->rErr += (r - t) + s;
   }
@@ -1701,11 +1701,11 @@
 }
 
 /*
 ** Add a (possibly large) integer to the running sum.
 */
-static void kahanBabuskaNeumaierStepInt64(volatile SumCtx *pSum, i64 iVal){
+static void kahanBabuskaNeumaierStepInt64(FP_VOLATILE SumCtx *pSum, i64 iVal){
   if( iVal<=-4503599627370496LL || iVal>=+4503599627370496LL ){
     i64 iBig, iSm;
     iSm = iVal % 16384;
     iBig = iVal - iSm;
     kahanBabuskaNeumaierStep(pSum, iBig);
@@ -1717,11 +1717,11 @@
 
 /*
 ** Initialize the Kahan-Babaska-Neumaier sum from a 64-bit integer
 */
 static void kahanBabuskaNeumaierInit(
-  volatile SumCtx *p,
+  FP_VOLATILE SumCtx *p,
   i64 iVal
 ){
   if( iVal<=-4503599627370496LL || iVal>=+4503599627370496LL ){
     i64 iSm = iVal % 16384;
     p->rSum = (double)(iVal - iSm);

Index: src/util.c
==================================================================
--- src/util.c
+++ src/util.c
@@ -390,20 +390,20 @@
 **
 ** Reference:
 **   T. J. Dekker, "A Floating-Point Technique for Extending the
 **   Available Precision".  1971-07-26.
 */
-static void dekkerMul2(volatile double *x, double y, double yy){
+static void dekkerMul2(FP_VOLATILE double *x, double y, double yy){
   /*
   ** The "volatile" keywords on parameter x[] and on local variables
   ** below are needed force intermediate results to be truncated to
   ** binary64 rather than be carried around in an extended-precision
   ** format.  The truncation is necessary for the Dekker algorithm to
   ** work.  Intel x86 floating point might omit the truncation without
   ** the use of volatile. 
   */
-  volatile double tx, ty, p, q, c, cc;
+  FP_VOLATILE double tx, ty, p, q, c, cc;
   double hx, hy;
   u64 m;
   memcpy(&m, (void*)&x[0], 8);
   m &= 0xfffffffffc000000LL;
   memcpy(&hx, &m, 8);

Why would that be "nice"? What does it accomplish? My gut is that it is better to avoid the added #ifdef complexity. But I'm open to opposing ideas, if you want to try to convince me.

As far s I know, that workaround is needed only for x86 and only if x87 fpu is used. If SSE is used for floating point calculations or code is compiled for any other architecture than x86 then the results are correct with or without the volatile.

The volatile prevents not only the incorrect optimizations but also the correct ones, unnecesarily degrading performance.

When compiling for 64-bit x86, compilers default to emit SSE instructions for floating point operations, and most compilers I've tested on godbolt.org even do it when compiling for 32-bit x86.

I just encountered a compiler that mishandles function parameters declared as volatile double, which affects one of the support routines for implementing the SUM/AVG/TOTAL aggregate functions (making them fail). Since the compiler vendor may view that usage of volatile as very rare, I don't know what priority they'll give to fixing their bug.

From what I can tell, the address of the variable in question never has its address taken so declaring it volatile is semantically meaningless in this situation.

Is there a unit test for checking the validity of the sumStep function and correct operation the Kahan-Babuska-Neumaier algorithm?

Here is a script you can feed into the CLI that tests sumStep using the Kahan-Babuska-Neumaier algorithm:

.echo on
.mode qbox
.testctrl uselongdouble off
CREATE TABLE t1(x REAL);
INSERT INTO t1(x) VALUES(1e50),(17.0),(-1e50);
SELECT sum(x) FROM t1;

That's useful. For the time being, I added an auto extension to override the built sum, total, and avg functions. It's a re-implementation with some re-factoring which avoid using volatile. Single letter variable names (e.g. r, s, t) may be the convention for academic work, but they are terrible for making understandable code.