@@ -2524,32 +2524,33 @@ With an argument, equivalent to object.__dict__.");
2524
2524
*/
2525
2525
2526
2526
typedef struct {
2527
- double sum ; /* accumulator */
2528
- double c; /* a running compensation for lost low-order bits */
2529
- } _csum ;
2527
+ double hi ; /* high-order bits for a running sum */
2528
+ double lo; /* a running compensation for lost low-order bits */
2529
+ } CompensatedSum ;
2530
2530
2531
2531
static inline void
2532
- _csum_neumaier_step(_csum *v, double x)
2532
+ cs_add(CompensatedSum *v, double x)
2533
2533
{
2534
- double t = v->sum + x;
2535
- if (fabs(v->sum ) >= fabs(x)) {
2536
- v->c += (v->sum - t) + x;
2534
+ double t = v->hi + x;
2535
+ if (fabs(v->hi ) >= fabs(x)) {
2536
+ v->lo += (v->hi - t) + x;
2537
2537
}
2538
2538
else {
2539
- v->c += (x - t) + v->sum ;
2539
+ v->lo += (x - t) + v->hi ;
2540
2540
}
2541
- v->sum = t;
2541
+ v->hi = t;
2542
2542
}
2543
2543
2544
- static inline void
2545
- _csum_neumaier_finalize(_csum *v)
2544
+ static inline double
2545
+ cs_to_double(CompensatedSum *v)
2546
2546
{
2547
2547
/* Avoid losing the sign on a negative result,
2548
2548
and don't let adding the compensation convert
2549
2549
an infinite or overflowed sum to a NaN. */
2550
- if (v->c && isfinite(v->c )) {
2551
- v->sum += v->c ;
2550
+ if (v->lo && isfinite(v->lo )) {
2551
+ v->hi += v->lo ;
2552
2552
}
2553
+ return v->hi;
2553
2554
}
2554
2555
2555
2556
/*[clinic input]
@@ -2664,19 +2665,18 @@ builtin_sum_impl(PyObject *module, PyObject *iterable, PyObject *start)
2664
2665
}
2665
2666
2666
2667
if (PyFloat_CheckExact(result)) {
2667
- _csum f_result = {PyFloat_AS_DOUBLE(result), 0.0 };
2668
+ CompensatedSum re_sum = {PyFloat_AS_DOUBLE(result)};
2668
2669
Py_SETREF(result, NULL);
2669
2670
while(result == NULL) {
2670
2671
item = PyIter_Next(iter);
2671
2672
if (item == NULL) {
2672
2673
Py_DECREF(iter);
2673
2674
if (PyErr_Occurred())
2674
2675
return NULL;
2675
- _csum_neumaier_finalize(&f_result);
2676
- return PyFloat_FromDouble(f_result.sum);
2676
+ return PyFloat_FromDouble(cs_to_double(&re_sum));
2677
2677
}
2678
2678
if (PyFloat_CheckExact(item)) {
2679
- _csum_neumaier_step(&f_result , PyFloat_AS_DOUBLE(item));
2679
+ cs_add(&re_sum , PyFloat_AS_DOUBLE(item));
2680
2680
_Py_DECREF_SPECIALIZED(item, _PyFloat_ExactDealloc);
2681
2681
continue;
2682
2682
}
@@ -2685,13 +2685,12 @@ builtin_sum_impl(PyObject *module, PyObject *iterable, PyObject *start)
2685
2685
int overflow;
2686
2686
value = PyLong_AsLongAndOverflow(item, &overflow);
2687
2687
if (!overflow) {
2688
- f_result.sum += (double)value;
2688
+ re_sum.hi += (double)value;
2689
2689
Py_DECREF(item);
2690
2690
continue;
2691
2691
}
2692
2692
}
2693
- _csum_neumaier_finalize(&f_result);
2694
- result = PyFloat_FromDouble(f_result.sum);
2693
+ result = PyFloat_FromDouble(cs_to_double(&re_sum));
2695
2694
if (result == NULL) {
2696
2695
Py_DECREF(item);
2697
2696
Py_DECREF(iter);
@@ -2710,8 +2709,8 @@ builtin_sum_impl(PyObject *module, PyObject *iterable, PyObject *start)
2710
2709
2711
2710
if (PyComplex_CheckExact(result)) {
2712
2711
Py_complex z = PyComplex_AsCComplex(result);
2713
- _csum cr_result = {z.real, 0.0 };
2714
- _csum ci_result = {z.imag, 0.0 };
2712
+ CompensatedSum re_sum = {z.real};
2713
+ CompensatedSum im_sum = {z.imag};
2715
2714
Py_SETREF(result, NULL);
2716
2715
while (result == NULL) {
2717
2716
item = PyIter_Next(iter);
@@ -2720,14 +2719,13 @@ builtin_sum_impl(PyObject *module, PyObject *iterable, PyObject *start)
2720
2719
if (PyErr_Occurred()) {
2721
2720
return NULL;
2722
2721
}
2723
- _csum_neumaier_finalize(&cr_result);
2724
- _csum_neumaier_finalize(&ci_result);
2725
- return PyComplex_FromDoubles(cr_result.sum, ci_result.sum);
2722
+ return PyComplex_FromDoubles(cs_to_double(&re_sum),
2723
+ cs_to_double(&im_sum));
2726
2724
}
2727
2725
if (PyComplex_CheckExact(item)) {
2728
2726
z = PyComplex_AsCComplex(item);
2729
- _csum_neumaier_step(&cr_result , z.real);
2730
- _csum_neumaier_step(&ci_result , z.imag);
2727
+ cs_add(&re_sum , z.real);
2728
+ cs_add(&im_sum , z.imag);
2731
2729
Py_DECREF(item);
2732
2730
continue;
2733
2731
}
@@ -2736,22 +2734,21 @@ builtin_sum_impl(PyObject *module, PyObject *iterable, PyObject *start)
2736
2734
int overflow;
2737
2735
value = PyLong_AsLongAndOverflow(item, &overflow);
2738
2736
if (!overflow) {
2739
- cr_result.sum += (double)value;
2740
- ci_result.sum += 0.0;
2737
+ re_sum.hi += (double)value;
2738
+ im_sum.hi += 0.0;
2741
2739
Py_DECREF(item);
2742
2740
continue;
2743
2741
}
2744
2742
}
2745
2743
if (PyFloat_Check(item)) {
2746
2744
double value = PyFloat_AS_DOUBLE(item);
2747
- cr_result.sum += value;
2748
- ci_result.sum += 0.0;
2745
+ re_sum.hi += value;
2746
+ im_sum.hi += 0.0;
2749
2747
_Py_DECREF_SPECIALIZED(item, _PyFloat_ExactDealloc);
2750
2748
continue;
2751
2749
}
2752
- _csum_neumaier_finalize(&cr_result);
2753
- _csum_neumaier_finalize(&ci_result);
2754
- result = PyComplex_FromDoubles(cr_result.sum, ci_result.sum);
2750
+ result = PyComplex_FromDoubles(cs_to_double(&re_sum),
2751
+ cs_to_double(&im_sum));
2755
2752
if (result == NULL) {
2756
2753
Py_DECREF(item);
2757
2754
Py_DECREF(iter);
0 commit comments