mirrors/sqlite
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Use the Kahan-Babushka-Neumaier algorithm to improve the accuracy of sum().

Browse Source 
FossilOrigin-Name: c63e26e705f5e967e14ef6aea8ce226548293ad8d25066069f29fa89673913d2
This commit is contained in:
drh 2023-07-06 15:44:38 +00:00
commit 59f13119cf
3 changed files with 80 additions and 22 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

14
manifest
View File

@ -1,5 +1,5 @@
C Use\s"volatile"\sisntead\sof\s"#pragma"\sto\sget\sfloating\spoint\scalculations\nworking\scorrectly\swhen\scompiling\swith\sGCC\sfor\sx86\smachines.
D 2023-07-06T00:55:06.161
C Use\sthe\sKahan-Babushka-Neumaier\salgorithm\sto\simprove\sthe\saccuracy\sof\ssum().
D 2023-07-06T15:44:38.951
F .fossil-settings/empty-dirs dbb81e8fc0401ac46a1491ab34a7f2c7c0452f2f06b54ebb845d024ca8283ef1
F .fossil-settings/ignore-glob 35175cdfcf539b2318cb04a9901442804be81cd677d8b889fcc9149c21f239ea
F LICENSE.md df5091916dbb40e6e9686186587125e1b2ff51f022cc334e886c19a0e9982724
@ -590,7 +590,7 @@ F src/delete.c cd5f5cd06ed0b6a882ec1a8c2a0d73b3cecb28479ad19e9931c4706c5e2182be
F src/expr.c 8d1656b65e26af3e34f78e947ac423f0d20c214ed25a67486e433bf16ca6b543
F src/fault.c 460f3e55994363812d9d60844b2a6de88826e007
F src/fkey.c a7fcbf7e66d14dbb73cf49f31489ebf66d0e6006c62b95246924a3bae9f37b36
F src/func.c 6028c160f693bdd018b651b5468a0a8e790f4e01e200796916b2d10a5d3237aa
F src/func.c 994a3a30d36c079c6e26d6653d770153675590daa549d3692eb38c5ac05dc663
F src/global.c a16553245e315ee0cda8f9b0bf744efef9dc99f86e9d77f58975ea58824ded92
F src/hash.c 9ee4269fb1d6632a6fecfb9479c93a1f29271bddbbaf215dd60420bcb80c7220
F src/hash.h 3340ab6e1d13e725571d7cee6d3e3135f0779a7d8e76a9ce0a85971fa3953c51
@ -2043,9 +2043,9 @@ F vsixtest/vsixtest.tcl 6a9a6ab600c25a91a7acc6293828957a386a8a93
F vsixtest/vsixtest.vcxproj.data 2ed517e100c66dc455b492e1a33350c1b20fbcdc
F vsixtest/vsixtest.vcxproj.filters 37e51ffedcdb064aad6ff33b6148725226cd608e
F vsixtest/vsixtest_TemporaryKey.pfx e5b1b036facdb453873e7084e1cae9102ccc67a0
P 7b4c16731e7bf6f03f5adf4fcb2008c0b19be473fb1b90b405c217c08916586a 1d972a690fdc70ab40862bd38427d68b48e8802ddf8e5c301f2d58ce2178b6ec
R edd6abc545b6b11f65503dd16b4a34c2
T +closed 1d972a690fdc70ab40862bd38427d68b48e8802ddf8e5c301f2d58ce2178b6ec
P 9427f42687ed6d97c474bf42d0c3e82d6f4b0075e74206adcb5699d72e32140e e3f7a960c9bd8e84cd70f0585bb955d043604a92001d0e2bf6c1216bb1fd7221
R fc07faedc8ec40abb250af878f96b89f
T +closed e3f7a960c9bd8e84cd70f0585bb955d043604a92001d0e2bf6c1216bb1fd7221
U drh
Z 9d851cc5730e074609430b016c224723
Z d5a624e7df2c6f31ea17c42fd16ebd48
# Remove this line to create a well-formed Fossil manifest.

2
manifest.uuid
View File

@ -1 +1 @@
9427f42687ed6d97c474bf42d0c3e82d6f4b0075e74206adcb5699d72e32140e
c63e26e705f5e967e14ef6aea8ce226548293ad8d25066069f29fa89673913d2

86
src/func.c
View File

@ -1671,12 +1671,58 @@ static void loadExt(sqlite3_context *context, int argc, sqlite3_value **argv){
typedef struct SumCtx SumCtx;
struct SumCtx {
double rSum; /* Running sum as as a double */
double rErr; /* Error term for Kahan-Babushka-Neumaier summation */
i64 iSum; /* Running sum as a signed integer */
i64 cnt; /* Number of elements summed */
u8 approx; /* True if any non-integer value was input to the sum */
u8 ovrfl; /* Integer overflow seen */
};
/*
** Do one step of the Kahan-Babushka-Neumaier summation.
**
** https://en.wikipedia.org/wiki/Kahan_summation_algorithm
**
** 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
){
volatile double s = pSum->rSum;
volatile double t = s + r;
if( fabs(s) > fabs(r) ){
pSum->rErr += (s - t) + r;
}else{
pSum->rErr += (r - t) + s;
}
pSum->rSum = t;
}
/*
** Add a (possibly large) integer to the running sum.
*/
static void kahanBabuskaNeumaierStepInt64(volatile SumCtx *pSum, i64 iVal){
volatile double rVal = (double)iVal;
kahanBabuskaNeumaierStep(pSum, rVal);
if( iVal<=-4503599627370496 || iVal>=+4503599627370496 ){
double rDiff = (double)(iVal - (i64)rVal);
kahanBabuskaNeumaierStep(pSum, rDiff);
}
}
/*
** Initialize the Kahan-Babaska-Neumaier sum from a 64-bit integer
*/
static void kahanBabuskaNeumaierInit(
volatile SumCtx *pSum,
i64 iVal
){
pSum->rSum = (double)iVal;
pSum->rErr = (double)(iVal - (i64)pSum->rSum);
}
/*
** Routines used to compute the sum, average, and total.
**
@ -1698,24 +1744,28 @@ static void sumStep(sqlite3_context *context, int argc, sqlite3_value **argv){
p->cnt++;
if( p->approx==0 ){
if( type!=SQLITE_INTEGER ){
p->rSum = (double)p->iSum;
kahanBabuskaNeumaierInit(p, p->iSum);
p->approx = 1;
p->rSum += sqlite3_value_double(argv[0]);
kahanBabuskaNeumaierStep(p, sqlite3_value_double(argv[0]));
}else{
i64 x = p->iSum;
if( sqlite3AddInt64(&x, sqlite3_value_int64(argv[0]))==0 ){
p->iSum = x;
}else{
p->ovrfl = 1;
p->rSum = (double)p->iSum;
kahanBabuskaNeumaierInit(p, p->iSum);
p->approx = 1;
p->rSum += sqlite3_value_double(argv[0]);
kahanBabuskaNeumaierStep(p, sqlite3_value_double(argv[0]));
}
}
}else{
if( type!=SQLITE_INTEGER ) p->ovrfl = 0;
p->approx = 1;
p->rSum += sqlite3_value_double(argv[0]);
if( type==SQLITE_INTEGER ){
kahanBabuskaNeumaierStepInt64(p, sqlite3_value_int64(argv[0]));
}else{
p->ovrfl = 0;
kahanBabuskaNeumaierStep(p, sqlite3_value_double(argv[0]));
}
}
}
}
@ -1732,10 +1782,18 @@ static void sumInverse(sqlite3_context *context, int argc, sqlite3_value**argv){
if( ALWAYS(p) && type!=SQLITE_NULL ){
assert( p->cnt>0 );
p->cnt--;
if( p->approx ){
p->rSum -= sqlite3_value_double(argv[0]);
}else{
if( !p->approx ){
p->iSum -= sqlite3_value_int64(argv[0]);
}else if( type==SQLITE_INTEGER ){
i64 iVal = sqlite3_value_int64(argv[0]);
if( iVal!=SMALLEST_INT64 ){
kahanBabuskaNeumaierStepInt64(p, -iVal);
}else{
kahanBabuskaNeumaierStepInt64(p, LARGEST_INT64);
kahanBabuskaNeumaierStepInt64(p, 1);
}
}else{
kahanBabuskaNeumaierStep(p, -sqlite3_value_double(argv[0]));
}
}
}
@ -1750,7 +1808,7 @@ static void sumFinalize(sqlite3_context *context){
if( p->ovrfl ){
sqlite3_result_error(context,"integer overflow",-1);
}else{
sqlite3_result_double(context, p->rSum);
sqlite3_result_double(context, p->rSum+p->rErr);
}
}else{
sqlite3_result_int64(context, p->iSum);
@ -1763,9 +1821,9 @@ static void avgFinalize(sqlite3_context *context){
if( p && p->cnt>0 ){
double r;
if( p->approx ){
r = p->rSum;
r = p->rSum+p->rErr;
}else{
r = sqlite3RealToI64(p->iSum);
r = (double)(p->iSum);
}
sqlite3_result_double(context, r/(double)p->cnt);
}
@ -1776,9 +1834,9 @@ static void totalFinalize(sqlite3_context *context){
p = sqlite3_aggregate_context(context, 0);
if( p ){
if( p->approx ){
r = p->rSum;
r = p->rSum+p->rErr;
}else{
r = sqlite3RealToI64(p->iSum);
r = (double)(p->iSum);
}
}
sqlite3_result_double(context, r);