| 1 | #ifndef DSTAR |
|---|
| 2 | |
|---|
| 3 | typedef unsigned int USItype; |
|---|
| 4 | typedef long int Wtype; |
|---|
| 5 | typedef long long int DWtype; |
|---|
| 6 | typedef long int word_type; |
|---|
| 7 | typedef unsigned long long int UDWtype; |
|---|
| 8 | typedef unsigned long int UWtype; |
|---|
| 9 | typedef unsigned int UQItype; |
|---|
| 10 | typedef double DFtype; |
|---|
| 11 | |
|---|
| 12 | #define BITS_PER_UNIT 8 |
|---|
| 13 | #define W_TYPE_SIZE (4 * BITS_PER_UNIT) |
|---|
| 14 | |
|---|
| 15 | #define __BITS4 (W_TYPE_SIZE / 4) |
|---|
| 16 | #define __ll_B ((UWtype) 1 << (W_TYPE_SIZE / 2)) |
|---|
| 17 | #define __ll_lowpart(t) ((UWtype) (t) & (__ll_B - 1)) |
|---|
| 18 | #define __ll_highpart(t) ((UWtype) (t) >> (W_TYPE_SIZE / 2)) |
|---|
| 19 | |
|---|
| 20 | #define udiv_qrnnd __udiv_qrnnd_c |
|---|
| 21 | #define __udiv_qrnnd_c(q, r, n1, n0, d) \ |
|---|
| 22 | do { \ |
|---|
| 23 | UWtype __d1, __d0, __q1, __q0; \ |
|---|
| 24 | UWtype __r1, __r0, __m; \ |
|---|
| 25 | __d1 = __ll_highpart (d); \ |
|---|
| 26 | __d0 = __ll_lowpart (d); \ |
|---|
| 27 | \ |
|---|
| 28 | __r1 = (n1) % __d1; \ |
|---|
| 29 | __q1 = (n1) / __d1; \ |
|---|
| 30 | __m = (UWtype) __q1 * __d0; \ |
|---|
| 31 | __r1 = __r1 * __ll_B | __ll_highpart (n0); \ |
|---|
| 32 | if (__r1 < __m) \ |
|---|
| 33 | { \ |
|---|
| 34 | __q1--, __r1 += (d); \ |
|---|
| 35 | if (__r1 >= (d)) /* i.e. we didn't get carry when adding to __r1 */\ |
|---|
| 36 | if (__r1 < __m) \ |
|---|
| 37 | __q1--, __r1 += (d); \ |
|---|
| 38 | } \ |
|---|
| 39 | __r1 -= __m; \ |
|---|
| 40 | \ |
|---|
| 41 | __r0 = __r1 % __d1; \ |
|---|
| 42 | __q0 = __r1 / __d1; \ |
|---|
| 43 | __m = (UWtype) __q0 * __d0; \ |
|---|
| 44 | __r0 = __r0 * __ll_B | __ll_lowpart (n0); \ |
|---|
| 45 | if (__r0 < __m) \ |
|---|
| 46 | { \ |
|---|
| 47 | __q0--, __r0 += (d); \ |
|---|
| 48 | if (__r0 >= (d)) \ |
|---|
| 49 | if (__r0 < __m) \ |
|---|
| 50 | __q0--, __r0 += (d); \ |
|---|
| 51 | } \ |
|---|
| 52 | __r0 -= __m; \ |
|---|
| 53 | \ |
|---|
| 54 | (q) = (UWtype) __q1 * __ll_B | __q0; \ |
|---|
| 55 | (r) = __r0; \ |
|---|
| 56 | } while (0) |
|---|
| 57 | |
|---|
| 58 | #define umul_ppmm(w1, w0, u, v) \ |
|---|
| 59 | __asm__ ("multu %2,%3" \ |
|---|
| 60 | : "=l" ((USItype) (w0)), \ |
|---|
| 61 | "=h" ((USItype) (w1)) \ |
|---|
| 62 | : "d" ((USItype) (u)), \ |
|---|
| 63 | "d" ((USItype) (v))) |
|---|
| 64 | #define UMUL_TIME 10 |
|---|
| 65 | #define UDIV_TIME 100 |
|---|
| 66 | |
|---|
| 67 | const UQItype __clz_tab[] = |
|---|
| 68 | { |
|---|
| 69 | 0,1,2,2,3,3,3,3,4,4,4,4,4,4,4,4,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5, |
|---|
| 70 | 6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6, |
|---|
| 71 | 7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7, |
|---|
| 72 | 7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7, |
|---|
| 73 | 8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8, |
|---|
| 74 | 8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8, |
|---|
| 75 | 8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8, |
|---|
| 76 | 8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8, |
|---|
| 77 | }; |
|---|
| 78 | |
|---|
| 79 | #define count_leading_zeros(count, x) \ |
|---|
| 80 | do { \ |
|---|
| 81 | UWtype __xr = (x); \ |
|---|
| 82 | UWtype __a; \ |
|---|
| 83 | \ |
|---|
| 84 | if (W_TYPE_SIZE <= 32) \ |
|---|
| 85 | { \ |
|---|
| 86 | __a = __xr < ((UWtype)1<<2*__BITS4) \ |
|---|
| 87 | ? (__xr < ((UWtype)1<<__BITS4) ? 0 : __BITS4) \ |
|---|
| 88 | : (__xr < ((UWtype)1<<3*__BITS4) ? 2*__BITS4 : 3*__BITS4); \ |
|---|
| 89 | } \ |
|---|
| 90 | else \ |
|---|
| 91 | { \ |
|---|
| 92 | for (__a = W_TYPE_SIZE - 8; __a > 0; __a -= 8) \ |
|---|
| 93 | if (((__xr >> __a) & 0xff) != 0) \ |
|---|
| 94 | break; \ |
|---|
| 95 | } \ |
|---|
| 96 | \ |
|---|
| 97 | (count) = W_TYPE_SIZE - (__clz_tab[__xr >> __a] + __a); \ |
|---|
| 98 | } while (0) |
|---|
| 99 | |
|---|
| 100 | struct DWstruct {Wtype low, high;}; |
|---|
| 101 | |
|---|
| 102 | typedef union |
|---|
| 103 | { |
|---|
| 104 | struct DWstruct s; |
|---|
| 105 | DWtype ll; |
|---|
| 106 | } DWunion; |
|---|
| 107 | |
|---|
| 108 | |
|---|
| 109 | #define sub_ddmmss(sh, sl, ah, al, bh, bl) \ |
|---|
| 110 | do { \ |
|---|
| 111 | UWtype __x; \ |
|---|
| 112 | __x = (al) - (bl); \ |
|---|
| 113 | (sh) = (ah) - (bh) - (__x > (al)); \ |
|---|
| 114 | (sl) = __x; \ |
|---|
| 115 | } while (0) |
|---|
| 116 | |
|---|
| 117 | |
|---|
| 118 | UDWtype __udivmoddi4 (UDWtype n, UDWtype d, UDWtype *rp) |
|---|
| 119 | { |
|---|
| 120 | const DWunion nn = {.ll = n}; |
|---|
| 121 | const DWunion dd = {.ll = d}; |
|---|
| 122 | DWunion rr; |
|---|
| 123 | UWtype d0, d1, n0, n1, n2; |
|---|
| 124 | UWtype q0, q1; |
|---|
| 125 | UWtype b, bm; |
|---|
| 126 | |
|---|
| 127 | d0 = dd.s.low; |
|---|
| 128 | d1 = dd.s.high; |
|---|
| 129 | n0 = nn.s.low; |
|---|
| 130 | n1 = nn.s.high; |
|---|
| 131 | |
|---|
| 132 | if (d1 == 0) |
|---|
| 133 | { |
|---|
| 134 | if (d0 > n1) |
|---|
| 135 | { |
|---|
| 136 | /* 0q = nn / 0D */ |
|---|
| 137 | |
|---|
| 138 | count_leading_zeros (bm, d0); |
|---|
| 139 | |
|---|
| 140 | if (bm != 0) |
|---|
| 141 | { |
|---|
| 142 | /* Normalize, i.e. make the most significant bit of the |
|---|
| 143 | denominator set. */ |
|---|
| 144 | |
|---|
| 145 | d0 = d0 << bm; |
|---|
| 146 | n1 = (n1 << bm) | (n0 >> (W_TYPE_SIZE - bm)); |
|---|
| 147 | n0 = n0 << bm; |
|---|
| 148 | } |
|---|
| 149 | |
|---|
| 150 | udiv_qrnnd (q0, n0, n1, n0, d0); |
|---|
| 151 | q1 = 0; |
|---|
| 152 | |
|---|
| 153 | /* Remainder in n0 >> bm. */ |
|---|
| 154 | } |
|---|
| 155 | else |
|---|
| 156 | { |
|---|
| 157 | /* qq = NN / 0d */ |
|---|
| 158 | |
|---|
| 159 | if (d0 == 0) |
|---|
| 160 | d0 = 1 / d0; /* Divide intentionally by zero. */ |
|---|
| 161 | |
|---|
| 162 | count_leading_zeros (bm, d0); |
|---|
| 163 | |
|---|
| 164 | if (bm == 0) |
|---|
| 165 | { |
|---|
| 166 | /* From (n1 >= d0) /\ (the most significant bit of d0 is set), |
|---|
| 167 | conclude (the most significant bit of n1 is set) /\ (the |
|---|
| 168 | leading quotient digit q1 = 1). |
|---|
| 169 | |
|---|
| 170 | This special case is necessary, not an optimization. |
|---|
| 171 | (Shifts counts of W_TYPE_SIZE are undefined.) */ |
|---|
| 172 | |
|---|
| 173 | n1 -= d0; |
|---|
| 174 | q1 = 1; |
|---|
| 175 | } |
|---|
| 176 | else |
|---|
| 177 | { |
|---|
| 178 | /* Normalize. */ |
|---|
| 179 | |
|---|
| 180 | b = W_TYPE_SIZE - bm; |
|---|
| 181 | |
|---|
| 182 | d0 = d0 << bm; |
|---|
| 183 | n2 = n1 >> b; |
|---|
| 184 | n1 = (n1 << bm) | (n0 >> b); |
|---|
| 185 | n0 = n0 << bm; |
|---|
| 186 | |
|---|
| 187 | udiv_qrnnd (q1, n1, n2, n1, d0); |
|---|
| 188 | } |
|---|
| 189 | |
|---|
| 190 | /* n1 != d0... */ |
|---|
| 191 | |
|---|
| 192 | udiv_qrnnd (q0, n0, n1, n0, d0); |
|---|
| 193 | |
|---|
| 194 | /* Remainder in n0 >> bm. */ |
|---|
| 195 | } |
|---|
| 196 | |
|---|
| 197 | if (rp != 0) |
|---|
| 198 | { |
|---|
| 199 | rr.s.low = n0 >> bm; |
|---|
| 200 | rr.s.high = 0; |
|---|
| 201 | *rp = rr.ll; |
|---|
| 202 | } |
|---|
| 203 | } |
|---|
| 204 | |
|---|
| 205 | else |
|---|
| 206 | { |
|---|
| 207 | if (d1 > n1) |
|---|
| 208 | { |
|---|
| 209 | /* 00 = nn / DD */ |
|---|
| 210 | |
|---|
| 211 | q0 = 0; |
|---|
| 212 | q1 = 0; |
|---|
| 213 | |
|---|
| 214 | /* Remainder in n1n0. */ |
|---|
| 215 | if (rp != 0) |
|---|
| 216 | { |
|---|
| 217 | rr.s.low = n0; |
|---|
| 218 | rr.s.high = n1; |
|---|
| 219 | *rp = rr.ll; |
|---|
| 220 | } |
|---|
| 221 | } |
|---|
| 222 | else |
|---|
| 223 | { |
|---|
| 224 | /* 0q = NN / dd */ |
|---|
| 225 | |
|---|
| 226 | count_leading_zeros (bm, d1); |
|---|
| 227 | if (bm == 0) |
|---|
| 228 | { |
|---|
| 229 | /* From (n1 >= d1) /\ (the most significant bit of d1 is set), |
|---|
| 230 | conclude (the most significant bit of n1 is set) /\ (the |
|---|
| 231 | quotient digit q0 = 0 or 1). |
|---|
| 232 | |
|---|
| 233 | This special case is necessary, not an optimization. */ |
|---|
| 234 | |
|---|
| 235 | /* The condition on the next line takes advantage of that |
|---|
| 236 | n1 >= d1 (true due to program flow). */ |
|---|
| 237 | if (n1 > d1 || n0 >= d0) |
|---|
| 238 | { |
|---|
| 239 | q0 = 1; |
|---|
| 240 | sub_ddmmss (n1, n0, n1, n0, d1, d0); |
|---|
| 241 | } |
|---|
| 242 | else |
|---|
| 243 | q0 = 0; |
|---|
| 244 | |
|---|
| 245 | q1 = 0; |
|---|
| 246 | |
|---|
| 247 | if (rp != 0) |
|---|
| 248 | { |
|---|
| 249 | rr.s.low = n0; |
|---|
| 250 | rr.s.high = n1; |
|---|
| 251 | *rp = rr.ll; |
|---|
| 252 | } |
|---|
| 253 | } |
|---|
| 254 | else |
|---|
| 255 | { |
|---|
| 256 | UWtype m1, m0; |
|---|
| 257 | /* Normalize. */ |
|---|
| 258 | |
|---|
| 259 | b = W_TYPE_SIZE - bm; |
|---|
| 260 | |
|---|
| 261 | d1 = (d1 << bm) | (d0 >> b); |
|---|
| 262 | d0 = d0 << bm; |
|---|
| 263 | n2 = n1 >> b; |
|---|
| 264 | n1 = (n1 << bm) | (n0 >> b); |
|---|
| 265 | n0 = n0 << bm; |
|---|
| 266 | |
|---|
| 267 | udiv_qrnnd (q0, n1, n2, n1, d1); |
|---|
| 268 | umul_ppmm (m1, m0, q0, d0); |
|---|
| 269 | |
|---|
| 270 | if (m1 > n1 || (m1 == n1 && m0 > n0)) |
|---|
| 271 | { |
|---|
| 272 | q0--; |
|---|
| 273 | sub_ddmmss (m1, m0, m1, m0, d1, d0); |
|---|
| 274 | } |
|---|
| 275 | |
|---|
| 276 | q1 = 0; |
|---|
| 277 | |
|---|
| 278 | /* Remainder in (n1n0 - m1m0) >> bm. */ |
|---|
| 279 | if (rp != 0) |
|---|
| 280 | { |
|---|
| 281 | sub_ddmmss (n1, n0, n1, n0, m1, m0); |
|---|
| 282 | rr.s.low = (n1 << b) | (n0 >> bm); |
|---|
| 283 | rr.s.high = n1 >> bm; |
|---|
| 284 | *rp = rr.ll; |
|---|
| 285 | } |
|---|
| 286 | } |
|---|
| 287 | } |
|---|
| 288 | } |
|---|
| 289 | |
|---|
| 290 | const DWunion ww = {{.low = q0, .high = q1}}; |
|---|
| 291 | return ww.ll; |
|---|
| 292 | } |
|---|
| 293 | |
|---|
| 294 | DWtype __divdi3 (DWtype u, DWtype v) |
|---|
| 295 | { |
|---|
| 296 | word_type c = 0; |
|---|
| 297 | DWunion uu = {.ll = u}; |
|---|
| 298 | DWunion vv = {.ll = v}; |
|---|
| 299 | DWtype w; |
|---|
| 300 | |
|---|
| 301 | if (uu.s.high < 0) |
|---|
| 302 | c = ~c, |
|---|
| 303 | uu.ll = -uu.ll; |
|---|
| 304 | if (vv.s.high < 0) |
|---|
| 305 | c = ~c, |
|---|
| 306 | vv.ll = -vv.ll; |
|---|
| 307 | |
|---|
| 308 | w = __udivmoddi4 (uu.ll, vv.ll, (UDWtype *) 0); |
|---|
| 309 | if (c) |
|---|
| 310 | w = -w; |
|---|
| 311 | |
|---|
| 312 | return w; |
|---|
| 313 | } |
|---|
| 314 | |
|---|
| 315 | DWtype __moddi3 (DWtype u, DWtype v) |
|---|
| 316 | { |
|---|
| 317 | word_type c = 0; |
|---|
| 318 | DWunion uu = {.ll = u}; |
|---|
| 319 | DWunion vv = {.ll = v}; |
|---|
| 320 | DWtype w; |
|---|
| 321 | |
|---|
| 322 | if (uu.s.high < 0) |
|---|
| 323 | c = ~c, |
|---|
| 324 | uu.ll = -uu.ll; |
|---|
| 325 | if (vv.s.high < 0) |
|---|
| 326 | vv.ll = -vv.ll; |
|---|
| 327 | |
|---|
| 328 | (void) __udivmoddi4 (uu.ll, vv.ll, &w); |
|---|
| 329 | if (c) |
|---|
| 330 | w = -w; |
|---|
| 331 | |
|---|
| 332 | return w; |
|---|
| 333 | } |
|---|
| 334 | |
|---|
| 335 | #define WORD_SIZE (sizeof (Wtype) * BITS_PER_UNIT) |
|---|
| 336 | #define HIGH_HALFWORD_COEFF (((UDWtype) 1) << (WORD_SIZE / 2)) |
|---|
| 337 | #define HIGH_WORD_COEFF (((UDWtype) 1) << WORD_SIZE) |
|---|
| 338 | |
|---|
| 339 | DFtype __floatdidf (DWtype u) |
|---|
| 340 | { |
|---|
| 341 | DFtype d = (Wtype) (u >> WORD_SIZE); |
|---|
| 342 | d *= HIGH_HALFWORD_COEFF; |
|---|
| 343 | d *= HIGH_HALFWORD_COEFF; |
|---|
| 344 | d += (UWtype) (u & (HIGH_WORD_COEFF - 1)); |
|---|
| 345 | |
|---|
| 346 | return d; |
|---|
| 347 | } |
|---|
| 348 | |
|---|
| 349 | DWtype __fixunsdfDI (DFtype a) |
|---|
| 350 | { |
|---|
| 351 | /* Get high part of result. The division here will just moves the radix |
|---|
| 352 | point and will not cause any rounding. Then the conversion to integral |
|---|
| 353 | type chops result as desired. */ |
|---|
| 354 | const UWtype hi = a / HIGH_WORD_COEFF; |
|---|
| 355 | |
|---|
| 356 | /* Get low part of result. Convert `hi' to floating type and scale it back, |
|---|
| 357 | then subtract this from the number being converted. This leaves the low |
|---|
| 358 | part. Convert that to integral type. */ |
|---|
| 359 | const UWtype lo = (a - ((DFtype) hi) * HIGH_WORD_COEFF); |
|---|
| 360 | |
|---|
| 361 | /* Assemble result from the two parts. */ |
|---|
| 362 | return ((UDWtype) hi << WORD_SIZE) | lo; |
|---|
| 363 | } |
|---|
| 364 | |
|---|
| 365 | DWtype __fixdfdi (DFtype a) |
|---|
| 366 | { |
|---|
| 367 | if (a < 0) |
|---|
| 368 | return - __fixunsdfDI (-a); |
|---|
| 369 | return __fixunsdfDI (a); |
|---|
| 370 | } |
|---|
| 371 | |
|---|
| 372 | UDWtype __umoddi3 (UDWtype u, UDWtype v) |
|---|
| 373 | { |
|---|
| 374 | UDWtype w; |
|---|
| 375 | |
|---|
| 376 | (void) __udivmoddi4 (u, v, &w); |
|---|
| 377 | |
|---|
| 378 | return w; |
|---|
| 379 | } |
|---|
| 380 | |
|---|
| 381 | UDWtype __udivdi3 (UDWtype n, UDWtype d) |
|---|
| 382 | { |
|---|
| 383 | return __udivmoddi4 (n, d, (UDWtype *) 0); |
|---|
| 384 | } |
|---|
| 385 | |
|---|
| 386 | #endif |
|---|