IsoSpec  2.1.2
summator.h
1 /*
2  * Copyright (C) 2015-2020 Mateusz Łącki and Michał Startek.
3  *
4  * This file is part of IsoSpec.
5  *
6  * IsoSpec is free software: you can redistribute it and/or modify
7  * it under the terms of the Simplified ("2-clause") BSD licence.
8  *
9  * IsoSpec is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
12  *
13  * You should have received a copy of the Simplified BSD Licence
14  * along with IsoSpec. If not, see <https://opensource.org/licenses/BSD-2-Clause>.
15  */
16 
17 #pragma once
18 
19 #include <cmath>
20 #include <vector>
21 #include <utility>
22 
23 namespace IsoSpec
24 {
25 
26 class SSummator
27 {
28  // Shewchuk algorithm
29  std::vector<double> partials;
30  int maxpart;
31  public:
32  inline SSummator()
33  { maxpart = 0; }
34 
35  inline SSummator(const SSummator& other) :
36  partials(other.partials),
37  maxpart(other.maxpart) {}
38 
39  inline void add(double x)
40  {
41  unsigned int i = 0;
42  for(int pidx = 0; pidx < maxpart; pidx++)
43  {
44  double y = partials[pidx];
45  if(std::abs(x) < std::abs(y))
46  std::swap(x, y);
47  double hi = x+y;
48  double lo = y-(hi-x);
49  if(lo != 0.0)
50  {
51  partials[i] = lo;
52  i += 1;
53  }
54  x = hi;
55  }
56  while(partials.size() <= i)
57  partials.push_back(0.0);
58  partials[i] = x;
59  maxpart = i+1;
60  }
61  inline double get()
62  {
63  double ret = 0.0;
64  for(int i = 0; i < maxpart; i++)
65  ret += partials[i];
66  return ret;
67  }
68 };
69 
70 
71 
72 
73 
74 
75 
76 class Summator{
77  // Kahan algorithm
78  double sum;
79  double c;
80 
81  public:
82  inline Summator()
83  { sum = 0.0; c = 0.0;}
84 
85  inline void add(double what)
86  {
87  double y = what - c;
88  double t = sum + y;
89  c = (t - sum) - y;
90  sum = t;
91  }
92 
93  inline double get()
94  {
95  return sum;
96  }
97 };
98 
99 class TSummator
100 {
101  // Trivial algorithm, for testing only
102  double sum;
103  public:
104  inline TSummator()
105  { sum = 0.0; }
106 
107  inline void add(double what)
108  {
109  sum += what;
110  }
111  inline double get()
112  {
113  return sum;
114  }
115 };
116 
117 
118 } // namespace IsoSpec
IsoSpec
Definition: allocator.cpp:20
IsoSpec::SSummator
Definition: summator.h:26
IsoSpec::TSummator
Definition: summator.h:99
IsoSpec::Summator
Definition: summator.h:76