diff options
Diffstat (limited to 'sum.c')
| -rw-r--r-- | sum.c | 26 |
1 files changed, 16 insertions, 10 deletions
@@ -1,6 +1,6 @@ /* Sum -- efficiently sum a list of floating-point numbers -Copyright 2004 Free Software Foundation, Inc. +Copyright 2004, 2005 Free Software Foundation, Inc. This file is part of the MPFR Library. @@ -67,7 +67,7 @@ void mpfr_count_sort (mpfr_ptr const tab[], unsigned long n, perm[target_rank] = tab[i]; } - TMP_FREE(account); + TMP_FREE(marker); } /* Performs a heap sort of the entries */ @@ -169,17 +169,18 @@ static int mpfr_list_sum_once (mpfr_ptr ret, mpfr_srcptr const tab[], { unsigned long i; mpfr_t sum; - int error_trap = 0; + int error_trap; + + if (MPFR_UNLIKELY (n == 1)) + return mpfr_set (ret, tab[0], GMP_RNDN); mpfr_init2 (sum, F); - mpfr_set (sum, tab[0], GMP_RNDN); + error_trap = mpfr_set (sum, tab[0], GMP_RNDN); for (i = 1; i < n - 1; i++) - { - error_trap |= mpfr_add (sum, sum, tab[i], GMP_RNDN); - } - + error_trap |= mpfr_add (sum, sum, tab[i], GMP_RNDN); error_trap |= mpfr_add (ret, sum, tab[n - 1], GMP_RNDN); + mpfr_clear (sum); return error_trap; } @@ -201,12 +202,17 @@ int mpfr_sum (mpfr_ptr ret, mpfr_ptr const tab[], unsigned long n, TMP_DECL(marker); TMP_MARK(marker); + if (MPFR_UNLIKELY (n == 0)) { + MPFR_SET_ZERO (ret); + MPFR_SET_POS (ret); + return 0; + } perm = (mpfr_srcptr *) TMP_ALLOC(n * sizeof(mpfr_srcptr)); mpfr_count_sort (tab, n, perm); - initial_f = MPFR_PREC(tab[0]); + initial_f = MAX (MPFR_PREC(tab[0]), MPFR_PREC(ret)); k = __gmpfr_ceil_log2 ((double) n) + 1; mpfr_init2 (cur_sum, initial_f); initial_guard_digits = k + 2; @@ -221,7 +227,7 @@ int mpfr_sum (mpfr_ptr ret, mpfr_ptr const tab[], unsigned long n, } while ((error_trap != 0) && !(mpfr_can_round (cur_sum, MPFR_GET_EXP(cur_sum) - current_f + 2, - GMP_RNDN, rnd, initial_f))); + GMP_RNDN, rnd, MPFR_PREC(ret)))); error_trap |= mpfr_set (ret, cur_sum, rnd); mpfr_clear (cur_sum); TMP_FREE(marker); |