Advertisement
homer512

Kahan

Jul 26th, 2014
294
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 1.15 KB | None | 0 0
  1. #include <iostream>
  2.  
  3. /*
  4.  * Demonstration of the Kahan summation algorithm
  5.  */
  6.  
  7. int main()
  8. {
  9.   float s = 0.;
  10.   float c = 0.;
  11.   for(long i = 1; i <= 100000000; ++i) {
  12.     float cur = 1.f / i / i;
  13.     float y = cur - c;
  14.     float t = s + y;
  15.     c = (t - s) - y;
  16.     s = t;
  17.   }
  18.   std::cout << s << std::endl;
  19.   s = c = 0.;
  20.   for(long i = 100000000; i >= 1; --i) {
  21.     float cur = 1.f / i / i;
  22.     float y = cur - c;
  23.     float t = s + y;
  24.     c = (t - s) - y;
  25.     s = t;
  26.   }
  27.   std::cout << s << std::endl;
  28. }
  29.  
  30. /*
  31.  * You can use this template to use it easier
  32.  */
  33.  
  34. template<class Real>
  35. class Accumulator
  36. {
  37.   Real _sum, _compensation;
  38. public:
  39.   Accumulator()
  40.     : _sum(0.), _compensation(0.)
  41.   {}
  42.   Real sum() const
  43.   { return _sum; }
  44.   Real operator()(Real added)
  45.   {
  46.     Real y = added - _compensation;
  47.     Real t = _sum + y;
  48.     _compensation = (t - _sum) - y;
  49.     _sum = t;
  50.     return _sum;
  51.   }
  52.   void reset(Real sum = 0.)
  53.   {
  54.     _sum = sum;
  55.     _compensation = 0.;
  56.   }
  57. };
  58.  
  59. void demonstration()
  60. {
  61.   Accumulator<float> add;
  62.   for(long i = 1; i <= 100000000; ++i)
  63.     add(1.f / i / i);
  64.   float sum = add.sum();
  65. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement