sum_kbn(fn, A) (issue 3, in part)

README.md

@@ -14,14 +14,14 @@ These functions were formerly part of Julia's Base library.
 ```julia
 julia> using KahanSummation
 
-julia> vec = [1.0, 1.0e16, -1.0e16, -0.5];
+julia> avec = [1.0, 1.0e16, -1.0e16, -0.5];
 
-julia> sum(vec), sum_kbn(vec)
+julia> sum(avec), sum_kbn(avec)
 (-0.5, 0.5)
 
-julia> vec = [1.0, 1.0e16, 1.0, -1.0e16];
+julia> avec = [1.0, 1.0e16, 1.0, -1.0e16];
 
-julia> cumsum_kbn(vec) .- cumsum(vec)
+julia> cumsum_kbn(avec) .- cumsum(avec)
 4-element Array{Float64,1}:
  0.0
  0.0

src/KahanSummation.jl

@@ -86,10 +86,11 @@ function cumsum_kbn(v::AbstractVector{T}) where T<:AbstractFloat
 end
 
 """
-    sum_kbn(A)
+    sum_kbn([f,] A)
 
 Return the sum of all elements of `A`, using the Kahan-Babuska-Neumaier compensated
-summation algorithm for additional accuracy.
+summation algorithm for additional accuracy.  When a function `f` is supplied, the
+result of calling f on each element of `A` is summed.
 """
 function sum_kbn(A)
     T = @default_eltype(typeof(A))
@@ -113,4 +114,26 @@ function sum_kbn(A)
     s - c
 end
 
+function sum_kbn(f,A)
+    T = @default_eltype(typeof(A))
+    c = promote_sys_size_add(zero(T)::T)
+    i = start(A)
+    if done(A, i)
+        return c
+    end
+    Ai, i = next(A, i)
+    s = f(Ai) - c
+    while !(done(A, i))
+        Ai, i = next(A, i)
+        fAi = f(Ai); t = s + fAi
+        if abs(s) >= abs(fAi)
+            c -= ((s-t) + fAi)
+        else
+            c -= ((fAi-t) + s)
+        end
+        s = t
+    end
+    s - c
+end
+
 end # module