This is the mail archive of the fortran@gcc.gnu.org mailing list for the GNU Fortran project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

Re: accuracy of intrinsic SUM function


On Fri, Aug 24, 2012 at 12:03:53PM +0200, Schorsch MCMLX wrote:
>
> The FNT95 error pattern exactly matches that of pairwise summation. 
> For a description of that algorithm, see e.g. wikipedia: 
> http://en.wikipedia.org/wiki/Pairwise_summation  . The gain in 
> accuracy resulting from using the pairwise algo instead of the 
> recurrent one is quite remarkable, as can be seen from the following 
> table (max. error in ULPs). 

(elided to keep this short) 

> Concluding, I'd like to propose to implement the pairwise summation 
> algorithm for the intrinsic SUM function in gFortran, if possible.


This looks quite interesting.  Do you timing measurements
and memory usage measurement?  At the moment, gfortran
uses the naive recurrence algorithm to accumulate the
sum.  gfortran in-line the summation, and thus creates
a single temporary variable to hold the sum.  Given

program test
   real x(10)
   x = [1., 2., 3., 4., 5., 6., 7., 8., 9., 0.]
   y = sum(x)
   print *, y
end program test

gfortran intermediate representation (ie., -fdump-tree-original)
gives

  real(kind=4) x[10];
  real(kind=4) y;

  {
    static real(kind=4) A.0[10] = {1.0e+0, 2.0e+0, 3.0e+0, 4.0e+0,\
          5.0e+0, 6.0e+0, 7.0e+0, 8.0e+0, 9.0e+0, 0.0};

    (void) (MEM[(c_char * {ref-all})&x] = MEM[(c_char * {ref-all})&A.0]);
    {
      real(kind=4) val.1;

      val.1 = 0.0;
      {
        integer(kind=8) S.2;

        S.2 = 1;
        while (1)
          {
            if (S.2 > 10) goto L.1;
            val.1 = x[S.2 + -1] + val.1;
            S.2 = S.2 + 1;
          }
        L.1:;
      }
      y = val.1;
    }
  
Will the recursive nature of the pairwise summation algorithm
generate memory pressure on the stack?  

Given the simplicist of the Kahan summation algorithm, it
may be better to replace gfortran naive implemenation 
with Kahan summation.  This would require 2 addditional
temporary variable and cause a slight increase execution
time. 

-- 
Steve


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]