☰
Go to the documentation of this file. 24 #ifndef _SCIP_DBLDBL_ARITH_ 25 #define _SCIP_DBLDBL_ARITH_ 30 #ifndef DISABLE_QUADPREC 33 #define QUAD_EPSILON 1e-12 36 #define QUAD_HI(x) x ## hi 37 #define QUAD_LO(x) x ## lo 38 #define QUAD(x) QUAD_HI(x), QUAD_LO(x) 39 #define QUAD_MEMBER(x) QUAD_HI(x); QUAD_LO(x) 40 #define QUAD_TO_DBL(x) ( QUAD_HI(x) + QUAD_LO(x) ) 41 #define QUAD_SCALE(x, a) do { QUAD_HI(x) *= (a); QUAD_LO(x) *= (a); } while(0) 42 #define QUAD_ASSIGN(a, constant) do { QUAD_HI(a) = (constant); QUAD_LO(a) = 0.0; } while(0) 43 #define QUAD_ASSIGN_Q(a, b) do { QUAD_HI(a) = QUAD_HI(b); QUAD_LO(a) = QUAD_LO(b); } while(0) 44 #define QUAD_ARRAY_SIZE(size) ((size)*2) 45 #define QUAD_ARRAY_LOAD(r, a, idx) do { QUAD_HI(r) = (a)[2*(idx)]; QUAD_LO(r) = (a)[2*(idx) + 1]; } while(0) 46 #define QUAD_ARRAY_STORE(a, idx, x) do { (a)[2*(idx)] = QUAD_HI(x); (a)[2*(idx) + 1] = QUAD_LO(x); } while(0) 49 #define SCIPquadprecProdDD(r, a, b) SCIPdbldblProd(QUAD_HI(r), QUAD_LO(r), a, b) 50 #define SCIPquadprecSquareD(r, a) SCIPdbldblSquare(QUAD_HI(r), QUAD_LO(r), a) 51 #define SCIPquadprecSumDD(r, a, b) SCIPdbldblSum(QUAD_HI(r), QUAD_LO(r), a, b) 52 #define SCIPquadprecDivDD(r, a, b) SCIPdbldblDiv(QUAD_HI(r), QUAD_LO(r), a, b) 53 #define SCIPquadprecSumQD(r, a, b) SCIPdbldblSum21(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a), b) 54 #define SCIPquadprecProdQD(r, a, b) SCIPdbldblProd21(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a), b) 55 #define SCIPquadprecDivDQ(r, a, b) SCIPdbldblDiv12(QUAD_HI(r), QUAD_LO(r), a, QUAD_HI(b), QUAD_LO(b)) 56 #define SCIPquadprecDivQD(r, a, b) SCIPdbldblDiv21(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a), b) 57 #define SCIPquadprecProdQQ(r, a, b) SCIPdbldblProd22(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a), QUAD_HI(b), QUAD_LO(b)) 58 #define SCIPquadprecSumQQ(r, a, b) SCIPdbldblSum22(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a), QUAD_HI(b), QUAD_LO(b)) 59 #define SCIPquadprecSquareQ(r, a) SCIPdbldblSquare2(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a)) 60 #define SCIPquadprecDivQQ(r, a, b) SCIPdbldblDiv22(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a), QUAD_HI(b), QUAD_LO(b)) 61 #define SCIPquadprecSqrtD(r, a) SCIPdbldblSqrt(QUAD_HI(r), QUAD_LO(r), a) 62 #define SCIPquadprecSqrtQ(r, a) SCIPdbldblSqrt2(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a)) 63 #define SCIPquadprecAbsQ(r, a) SCIPdbldblAbs2(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a)) 64 #define SCIPquadprecFloorQ(r, a) SCIPdbldblFloor2(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a)) 65 #define SCIPquadprecCeilQ(r, a) SCIPdbldblCeil2(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a)) 66 #define SCIPquadprecEpsFloorQ(r, a, eps) SCIPdbldblEpsFloor2(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a), eps) 67 #define SCIPquadprecEpsCeilQ(r, a, eps) SCIPdbldblEpsCeil2(QUAD_HI(r), QUAD_LO(r), QUAD_HI(a), QUAD_LO(a), eps) 72 #define QUAD_EPSILON 1e-9 76 #define QUAD_LO(x) 0.0 78 #define QUAD_MEMBER(x) x 79 #define QUAD_TO_DBL(x) (x) 80 #define QUAD_SCALE(x, a) do { (x) *= (a); } while(0) 81 #define QUAD_ASSIGN(a, constant) do { (a) = constant; } while(0) 82 #define QUAD_ASSIGN_Q(a, b) do { (a) = (b); } while(0) 83 #define QUAD_ARRAY_SIZE(size) (size) 84 #define QUAD_ARRAY_LOAD(r, a, idx) do { r = (a)[(idx)]; } while(0) 85 #define QUAD_ARRAY_STORE(a, idx, x) do { (a)[(idx)] = (x); } while(0) 87 #define SCIPquadprecProdDD(r, a, b) do { (r) = (a) * (b); } while(0) 88 #define SCIPquadprecSquareD(r, a) do { (r) = (a) * (a); } while(0) 89 #define SCIPquadprecSumDD(r, a, b) do { (r) = (a) + (b); } while(0) 90 #define SCIPquadprecDivDD(r, a, b) do { (r) = (a) / (b); } while(0) 91 #define SCIPquadprecSumQD(r, a, b) do { (r) = (a) + (b); } while(0) 92 #define SCIPquadprecProdQD(r, a, b) do { (r) = (a) * (b); } while(0) 93 #define SCIPquadprecDivDQ(r, a, b) do { (r) = (a) / (b); } while(0) 94 #define SCIPquadprecDivQD(r, a, b) do { (r) = (a) / (b); } while(0) 95 #define SCIPquadprecProdQQ(r, a, b) do { (r) = (a) * (b); } while(0) 96 #define SCIPquadprecSumQQ(r, a, b) do { (r) = (a) + (b); } while(0) 97 #define SCIPquadprecSquareQ(r, a) do { (r) = (a) * (a); } while(0) 98 #define SCIPquadprecDivQQ(r, a, b) do { (r) = (a) / (b); } while(0) 99 #define SCIPquadprecSqrtD(r, a) do { (r) = sqrt(a); } while(0) 100 #define SCIPquadprecSqrtQ(r, a) do { (r) = sqrt(a); } while(0) 101 #define SCIPquadprecAbsQ(r, a) do { (r) = fabs(a); } while(0) 102 #define SCIPquadprecFloorQ(r, a) do { (r) = floor(a); } while(0) 103 #define SCIPquadprecCeilQ(r, a) do { (r) = ceil(a); } while(0) 104 #define SCIPquadprecEpsFloorQ(r, a, eps) do { (r) = floor((a) + (eps)); } while(0) 105 #define SCIPquadprecEpsCeilQ(r, a, eps) do { (r) = ceil((a) - (eps)); } while(0) 109 #define __SCIPdbldblSplit(rhi, rlo, x) \ 111 const double __tmp_split_dbl = 134217729.0 * (x); \ 112 (rhi) = __tmp_split_dbl - (__tmp_split_dbl - (x)); \ 113 (rlo) = (x) - (rhi);\ 117 #define SCIPdbldblProd(rhi, rlo, a, b) \ 119 double __tmp_dbldbl_prod_ahi; \ 120 double __tmp_dbldbl_prod_alo; \ 121 double __tmp_dbldbl_prod_bhi; \ 122 double __tmp_dbldbl_prod_blo; \ 123 __SCIPdbldblSplit(__tmp_dbldbl_prod_ahi, __tmp_dbldbl_prod_alo, a); \ 124 __SCIPdbldblSplit(__tmp_dbldbl_prod_bhi, __tmp_dbldbl_prod_blo, b); \ 126 (rlo) = __tmp_dbldbl_prod_alo * __tmp_dbldbl_prod_blo - \ 127 ((((rhi) - __tmp_dbldbl_prod_ahi * __tmp_dbldbl_prod_bhi) \ 128 - __tmp_dbldbl_prod_alo * __tmp_dbldbl_prod_bhi) \ 129 - __tmp_dbldbl_prod_ahi * __tmp_dbldbl_prod_blo); \ 133 #define SCIPdbldblSquare(rhi, rlo, a) \ 135 double __tmp_dbldbl_square_ahi; \ 136 double __tmp_dbldbl_square_alo; \ 137 __SCIPdbldblSplit(__tmp_dbldbl_square_ahi, __tmp_dbldbl_square_alo, a); \ 139 (rlo) = __tmp_dbldbl_square_alo * __tmp_dbldbl_square_alo - \ 140 ((((rhi) - __tmp_dbldbl_square_ahi * __tmp_dbldbl_square_ahi) \ 141 - 2.0 * __tmp_dbldbl_square_alo * __tmp_dbldbl_square_ahi)); \ 145 #define SCIPdbldblSum(rhi, rlo, a, b) \ 147 double __tmp1_dbldbl_sum; \ 148 double __tmp2_dbldbl_sum; \ 149 __tmp2_dbldbl_sum = (a) + (b); \ 150 __tmp1_dbldbl_sum = __tmp2_dbldbl_sum - (a); \ 151 (rlo) = ((a) - (__tmp2_dbldbl_sum - __tmp1_dbldbl_sum)) + ((b) - __tmp1_dbldbl_sum); \ 152 (rhi) = __tmp2_dbldbl_sum; \ 156 #define SCIPdbldblDiv(rhi, rlo, a, b) \ 158 double __tmp_dbldbl_div_hi; \ 159 double __tmp_dbldbl_div_lo; \ 160 double __estim_dbldbl_div = (a)/(b); \ 161 SCIPdbldblProd(__tmp_dbldbl_div_hi, __tmp_dbldbl_div_lo, b, __estim_dbldbl_div); \ 162 SCIPdbldblSum21(__tmp_dbldbl_div_hi, __tmp_dbldbl_div_lo, __tmp_dbldbl_div_hi, __tmp_dbldbl_div_lo, -(a)); \ 163 __tmp_dbldbl_div_hi /= (b); \ 164 __tmp_dbldbl_div_lo /= (b); \ 165 SCIPdbldblSum21(rhi, rlo, -__tmp_dbldbl_div_hi, -__tmp_dbldbl_div_lo, __estim_dbldbl_div); \ 171 #define SCIPdbldblSum21(rhi, rlo, ahi, alo, b) \ 173 double __tmp_dbldbl_sum21_hi; \ 174 double __tmp_dbldbl_sum21_lo; \ 175 SCIPdbldblSum(__tmp_dbldbl_sum21_hi, __tmp_dbldbl_sum21_lo, ahi, b); \ 176 (rlo) = __tmp_dbldbl_sum21_lo + (alo); \ 177 (rhi) = __tmp_dbldbl_sum21_hi; \ 184 #define SCIPdbldblProd21(rhi, rlo, ahi, alo, b) \ 186 double __tmp_dbldbl_prod21_hi; \ 187 double __tmp_dbldbl_prod21_lo; \ 188 SCIPdbldblProd(__tmp_dbldbl_prod21_hi, __tmp_dbldbl_prod21_lo, ahi, b); \ 189 (rlo) = (alo) * (b) + __tmp_dbldbl_prod21_lo; \ 190 (rhi) = __tmp_dbldbl_prod21_hi; \ 196 #define SCIPdbldblDiv12(rhi, rlo, a, bhi, blo) \ 198 double __tmp_dbldbl_div12_hi; \ 199 double __tmp_dbldbl_div12_lo; \ 200 double __estim_dbldbl_div12 = (a)/(bhi); \ 201 SCIPdbldblProd21(__tmp_dbldbl_div12_hi, __tmp_dbldbl_div12_lo, bhi, blo, __estim_dbldbl_div12); \ 202 SCIPdbldblSum21(__tmp_dbldbl_div12_hi, __tmp_dbldbl_div12_lo, __tmp_dbldbl_div12_hi, __tmp_dbldbl_div12_lo, -(a)); \ 203 __tmp_dbldbl_div12_hi /= (bhi); \ 204 __tmp_dbldbl_div12_lo /= (bhi); \ 205 SCIPdbldblSum21(rhi, rlo, -__tmp_dbldbl_div12_hi, -__tmp_dbldbl_div12_lo, __estim_dbldbl_div12); \ 212 #define SCIPdbldblDiv21(rhi, rlo, ahi, alo, b) \ 214 double __tmp_dbldbl_div21_hi; \ 215 double __tmp_dbldbl_div21_lo; \ 216 double __estim_dbldbl_div21_hi; \ 217 double __estim_dbldbl_div21_lo; \ 218 __estim_dbldbl_div21_hi = (ahi)/(b); \ 219 __estim_dbldbl_div21_lo = (alo)/(b); \ 220 SCIPdbldblProd21(__tmp_dbldbl_div21_hi, __tmp_dbldbl_div21_lo, __estim_dbldbl_div21_hi, __estim_dbldbl_div21_lo, b); \ 221 SCIPdbldblSum22(__tmp_dbldbl_div21_hi, __tmp_dbldbl_div21_lo, __tmp_dbldbl_div21_hi, __tmp_dbldbl_div21_lo, -(ahi), -(alo)); \ 222 __tmp_dbldbl_div21_hi /= (b); \ 223 __tmp_dbldbl_div21_lo /= (b); \ 224 SCIPdbldblSum22(rhi, rlo, __estim_dbldbl_div21_hi, __estim_dbldbl_div21_lo, -__tmp_dbldbl_div21_hi, -__tmp_dbldbl_div21_lo); \ 228 #define SCIPdbldblProd22(rhi, rlo, ahi, alo, bhi, blo) \ 230 double __tmp_dbldbl_prod22_hi; \ 231 double __tmp_dbldbl_prod22_lo; \ 232 SCIPdbldblProd(__tmp_dbldbl_prod22_hi, __tmp_dbldbl_prod22_lo, ahi, bhi); \ 233 SCIPdbldblSum21(__tmp_dbldbl_prod22_hi, __tmp_dbldbl_prod22_lo, \ 234 __tmp_dbldbl_prod22_hi, __tmp_dbldbl_prod22_lo, (alo) * (bhi)); \ 235 SCIPdbldblSum21(rhi, rlo, \ 236 __tmp_dbldbl_prod22_hi, __tmp_dbldbl_prod22_lo, (ahi) * (blo)); \ 240 #define SCIPdbldblSum22(rhi, rlo, ahi, alo, bhi, blo) \ 242 double __tmp_dbldbl_sum22_hi; \ 243 double __tmp_dbldbl_sum22_lo; \ 244 SCIPdbldblSum21(__tmp_dbldbl_sum22_hi, __tmp_dbldbl_sum22_lo, ahi, alo, bhi); \ 245 SCIPdbldblSum21(rhi, rlo, __tmp_dbldbl_sum22_hi, __tmp_dbldbl_sum22_lo, blo); \ 249 #define SCIPdbldblSquare2(rhi, rlo, ahi, alo) \ 251 double __tmp_dbldbl_square2_hi; \ 252 double __tmp_dbldbl_square2_lo; \ 253 SCIPdbldblSquare(__tmp_dbldbl_square2_hi, __tmp_dbldbl_square2_lo, (ahi)); \ 254 SCIPdbldblSum21(rhi, rlo, __tmp_dbldbl_square2_hi, __tmp_dbldbl_square2_lo, 2 * (ahi) * (alo)); \ 258 #define SCIPdbldblDiv22(rhi, rlo, ahi, alo, bhi, blo) \ 260 double __tmp_dbldbl_div22_hi; \ 261 double __tmp_dbldbl_div22_lo; \ 262 double __estim_dbldbl_div22_hi = (ahi) / (bhi); \ 263 double __estim_dbldbl_div22_lo = (alo) / (bhi); \ 264 SCIPdbldblProd22(__tmp_dbldbl_div22_hi, __tmp_dbldbl_div22_lo, \ 265 bhi, blo, __estim_dbldbl_div22_hi, __estim_dbldbl_div22_lo); \ 266 SCIPdbldblSum22(__tmp_dbldbl_div22_hi, __tmp_dbldbl_div22_lo, \ 267 __tmp_dbldbl_div22_hi, __tmp_dbldbl_div22_lo, -(ahi), -(alo)); \ 268 __tmp_dbldbl_div22_hi /= (bhi); \ 269 __tmp_dbldbl_div22_lo /= (bhi); \ 270 SCIPdbldblSum22(rhi, rlo, __estim_dbldbl_div22_hi, __estim_dbldbl_div22_lo, \ 271 -__tmp_dbldbl_div22_hi, -__tmp_dbldbl_div22_lo); \ 276 #define SCIPdbldblSqrt(rhi, rlo, a) \ 278 double __estim_dbldbl_sqrt = sqrt(a); \ 279 if( __estim_dbldbl_sqrt != 0.0 ) \ 281 SCIPdbldblDiv(rhi, rlo, a, __estim_dbldbl_sqrt); \ 282 SCIPdbldblSum21(rhi, rlo, rhi, rlo, __estim_dbldbl_sqrt); \ 295 #define SCIPdbldblSqrt2(rhi, rlo, ahi, alo) \ 297 double __estim_dbldbl_sqrt2 = sqrt(ahi + alo); \ 298 if( __estim_dbldbl_sqrt2 != 0.0 ) \ 300 SCIPdbldblDiv21(rhi, rlo, ahi, alo, __estim_dbldbl_sqrt2); \ 301 SCIPdbldblSum21(rhi, rlo, rhi, rlo, __estim_dbldbl_sqrt2); \ 313 #define SCIPdbldblAbs2(rhi, rlo, ahi, alo) \ 328 #define SCIPdbldblFloor2(rhi, rlo, ahi, alo) \ 330 double __tmp_dbldbl_floor; \ 331 __tmp_dbldbl_floor = floor((ahi) + (alo)); \ 332 SCIPdbldblSum21(rhi, rlo, ahi, alo, -__tmp_dbldbl_floor); \ 333 if( ((rhi) - 1.0) + (rlo) < 0.0 && (rhi) + (rlo) >= 0.0 ) \ 336 (rhi) = __tmp_dbldbl_floor; \ 342 double __tmp2_dbldbl_floor = floor((rhi) + (rlo)); \ 343 SCIPdbldblSum(rhi, rlo, __tmp_dbldbl_floor, __tmp2_dbldbl_floor); \ 348 #define SCIPdbldblCeil2(rhi, rlo, ahi, alo) \ 350 double __tmp_dbldbl_ceil; \ 351 __tmp_dbldbl_ceil = ceil((ahi) + (alo)); \ 352 SCIPdbldblSum21(rhi, rlo, -(ahi), -(alo), __tmp_dbldbl_ceil); \ 353 if( ((rhi) - 1.0) + (rlo) < 0.0 && (rhi) + (rlo) >= 0.0 ) \ 356 (rhi) = __tmp_dbldbl_ceil; \ 362 double __tmp2_dbldbl_ceil = floor((rhi) + (rlo)); \ 363 SCIPdbldblSum(rhi, rlo, __tmp_dbldbl_ceil, -__tmp2_dbldbl_ceil); \ 368 #define SCIPdbldblEpsFloor2(rhi, rlo, ahi, alo, eps) SCIPdbldblFloor2(rhi, rlo, ahi, (alo) + (eps)) 371 #define SCIPdbldblEpsCeil2(rhi, rlo, ahi, alo, eps) SCIPdbldblCeil2(rhi, rlo, ahi, (alo) - (eps))