|
45 | 45 |
|
46 | 46 | #include "pgs_util.h"
|
47 | 47 |
|
| 48 | +/* On some 32 bit platforms, there is a gcc bug that makes floating point |
| 49 | + * calculations and comparisons unstable (see the link below). The problem |
| 50 | + * originates in FPU 80 bits registers where double values are not truncated |
| 51 | + * to 64 bit values. When gcc compiles some code with enabled optimizations, |
| 52 | + * the intermediate results may be kept in the FPU registers without truncation |
| 53 | + * to 64 bit values. Extra bits may produce unstable results when comparing |
| 54 | + * the numbers. |
| 55 | + * |
| 56 | + * The generic solution is to save the intermediate results in the memory where |
| 57 | + * the values are truncated to 64 bit values. It affects the performance but |
| 58 | + * makes the tests stable on all platforms. |
| 59 | + * |
| 60 | + * PGSPHERE_FLOAT_STORE macro enables storing of intermediate results for FPxx |
| 61 | + * operations in the memory. It is enabled by default for 32 bit platforms. |
| 62 | + * It can be explicitly enabled or disabled in CFLAGS. To enable it for all |
| 63 | + * code the gcc option -ffloat-store may be used as well. |
| 64 | + * |
| 65 | + * Link to gcc bug: https://gcc.gnu.org/bugzilla/show_bug.cgi?id=323 |
| 66 | + */ |
| 67 | +#if !defined(PGSPHERE_FLOAT_STORE) |
| 68 | +#if _WIN64 || (__GNUC__ && __x86_64__) |
| 69 | +#define PGSPHERE_FLOAT_STORE 0 |
| 70 | +#elif _WIN32 || __GNUC__ |
| 71 | +#define PGSPHERE_FLOAT_STORE 1 |
| 72 | +#else |
| 73 | +#define PGSPHERE_FLOAT_STORE 0 |
| 74 | +#endif |
| 75 | +#endif // PGSPHERE_FLOAT_STORE |
| 76 | + |
48 | 77 | #define EPSILON 1.0E-09
|
49 | 78 |
|
50 | 79 | #define FPzero(A) (fabs(A) <= EPSILON)
|
51 | 80 |
|
| 81 | +#if PGSPHERE_FLOAT_STORE |
| 82 | + |
| 83 | +static inline bool |
| 84 | +FPeq(double A, double B) |
| 85 | +{ |
| 86 | + const volatile double AB = A - B; |
| 87 | + return A == B || fabs(AB) <= EPSILON; |
| 88 | +} |
| 89 | + |
| 90 | +static inline bool |
| 91 | +FPne(double A, double B) |
| 92 | +{ |
| 93 | + const volatile double AB = A - B; |
| 94 | + return A != B && fabs(AB) > EPSILON; |
| 95 | +} |
| 96 | + |
| 97 | +static inline bool |
| 98 | +FPlt(double A, double B) |
| 99 | +{ |
| 100 | + const volatile double AE = A + EPSILON; |
| 101 | + return AE < B; |
| 102 | +} |
| 103 | + |
| 104 | +static inline bool |
| 105 | +FPle(double A, double B) |
| 106 | +{ |
| 107 | + const volatile double BE = B + EPSILON; |
| 108 | + return A <= BE; |
| 109 | +} |
| 110 | + |
| 111 | +static inline bool |
| 112 | +FPgt(double A, double B) |
| 113 | +{ |
| 114 | + const volatile double BE = B + EPSILON; |
| 115 | + return A > BE; |
| 116 | +} |
| 117 | + |
| 118 | +static inline bool |
| 119 | +FPge(double A, double B) |
| 120 | +{ |
| 121 | + const volatile double AE = A + EPSILON; |
| 122 | + return AE >= B; |
| 123 | +} |
| 124 | + |
| 125 | +#else |
| 126 | + |
52 | 127 | static inline bool
|
53 | 128 | FPeq(double A, double B)
|
54 | 129 | {
|
@@ -85,6 +160,8 @@ FPge(double A, double B)
|
85 | 160 | return A + EPSILON >= B;
|
86 | 161 | }
|
87 | 162 |
|
| 163 | +#endif // PGSPHERE_FLOAT_STORE |
| 164 | + |
88 | 165 | /*---------------------------------------------------------------------
|
89 | 166 | * Point - (x,y)
|
90 | 167 | *-------------------------------------------------------------------*/
|
|
0 commit comments