HiRep 0.1
Loading...
Searching...
No Matches
precise_sums.h
1/***************************************************************************\
2* Copyright (c) 2022, Claudio Pica *
3* All rights reserved. *
4\***************************************************************************/
5
6#ifndef PRECISE_SUMS_H
7#define PRECISE_SUMS_H
8
9#include <math.h>
10
12#define _comp_term(s, a, b) (fabs(a) > fabs(b)) ? (((a) - (s)) + (b)) : (((b) - (s)) + (a))
13#define _2Sum_acc(s, c, a) _2Sum_acc_hirep(s, c, a)
14
18#define _2Sum_acc_hirep(s, c, a) \
19 do { \
20 double _a = a; \
21 double _t = s + _a; \
22 c += _comp_term(_t, s, _a); \
23 s = _t; \
24 } while (0)
25
27#define _2Sum_acc_hirep2(s, c, cc, a) \
28 do { \
29 double _a = a; \
30 double _t = s + _a; \
31 double _ct = _comp_term(_t, s, _a); \
32 s = _t; \
33 _t = c + _ct; \
34 double _cct = _comp_term(_t, c, _ct); \
35 c = _t; \
36 cc += _cct; \
37 } while (0)
38
40#define _2Sum_acc_oqcd(s, c, a) \
41 do { \
42 double _a = s + a; \
43 double _qp = _a - a; \
44 double _up = _a - _qp; \
45 double _b = (s - _qp) + (a - _up); \
46 double _c = c + _b; \
47 double _d = _a + _c; \
48 s = _d; \
49 c = _c - (_d - _a); \
50 } while (0)
51
52#define _2Sum_acc_naive(s, c, a) \
53 do { \
54 (s) += (a); \
55 c = 0.; \
56 } while (0)
57
58#endif