comparison src/number-mp.c @ 5602:c9e5612f5424

Support the MP library on recent FreeBSD, have it pass relevant tests. src/ChangeLog addition: 2011-11-26 Aidan Kehoe <kehoea@parhasard.net> * number-mp.c (bignum_to_string): Don't overwrite the accumulator we've just set up for this function. * number-mp.c (BIGNUM_TO_TYPE): mp_itom() doesn't necessarily do what this code used to think with negative numbers, it can treat them as unsigned ints. Subtract numbers from bignum_zero instead of multiplying them by -1 to convert them to their negative equivalents. * number-mp.c (bignum_to_int): * number-mp.c (bignum_to_uint): * number-mp.c (bignum_to_long): * number-mp.c (bignum_to_ulong): * number-mp.c (bignum_to_double): Use the changed BIGNUM_TO_TYPE() in these functions. * number-mp.c (bignum_ceil): * number-mp.c (bignum_floor): In these functions, be more careful about rounding to positive and negative infinity, respectively. Don't use the sign of QUOTIENT when working out out whether to add or subtract one, rather use the sign QUOTIENT would have if arbitrary-precision division were done. * number-mp.h: * number-mp.h (MP_GCD): Wrap #include <mp.h> in BEGIN_C_DECLS/END_C_DECLS. * number.c (Fbigfloat_get_precision): * number.c (Fbigfloat_set_precision): Don't attempt to call XBIGFLOAT_GET_PREC if this build doesn't support big floats.
author Aidan Kehoe <kehoea@parhasard.net>
date Sat, 26 Nov 2011 17:59:14 +0000
parents 2aa9cd456ae7
children 7aa144d1404b
comparison
equal deleted inserted replaced
5601:3e5d5e8e4bb7 5602:c9e5612f5424
40 short rem; 40 short rem;
41 41
42 /* FIXME: signal something if base is < 2 or doesn't fit into a short. */ 42 /* FIXME: signal something if base is < 2 or doesn't fit into a short. */
43 43
44 /* Save the sign for later */ 44 /* Save the sign for later */
45 sign = MP_MCMP (b, bignum_zero); 45 sign = bignum_sign (b);
46 46
47 if (sign == 0) 47 if (sign == 0)
48 { 48 {
49 XREALLOC_ARRAY (buffer, char, 2); 49 XREALLOC_ARRAY (buffer, char, 2);
50 buffer[0] = '0'; 50 buffer[0] = '0';
54 /* Copy abs(b) into quo for destructive modification */ 54 /* Copy abs(b) into quo for destructive modification */
55 else if (sign < 0) 55 else if (sign < 0)
56 MP_MSUB (bignum_zero, b, quo); 56 MP_MSUB (bignum_zero, b, quo);
57 else 57 else
58 MP_MOVE (b, quo); 58 MP_MOVE (b, quo);
59
60 quo = MP_ITOM (0);
61 59
62 /* Loop over the digits of b (in BASE) and place each one into buffer */ 60 /* Loop over the digits of b (in BASE) and place each one into buffer */
63 for (i = 0U; MP_MCMP(quo, bignum_zero) > 0; i++) 61 for (i = 0U; MP_MCMP(quo, bignum_zero) > 0; i++)
64 { 62 {
65 MP_SDIV (quo, base, quo, &rem); 63 MP_SDIV (quo, base, quo, &rem);
88 xfree (buffer); 86 xfree (buffer);
89 return retval; 87 return retval;
90 } 88 }
91 89
92 #define BIGNUM_TO_TYPE(type,accumtype) do { \ 90 #define BIGNUM_TO_TYPE(type,accumtype) do { \
93 MP_MULT (b, quo, quo); \ 91 if (0 == sign) \
92 { \
93 return (type)0; \
94 } \
95 \
96 bignum_init (quo); \
97 \
98 if (sign < 0) \
99 { \
100 MP_MSUB (bignum_zero, b, quo); \
101 } \
102 else \
103 { \
104 MP_MOVE (b, quo); \
105 } \
106 \
94 for (i = 0U; i < sizeof(type); i++) \ 107 for (i = 0U; i < sizeof(type); i++) \
95 { \ 108 { \
96 MP_SDIV (quo, 256, quo, &rem); \ 109 MP_SDIV (quo, 256, quo, &rem); \
97 retval |= ((accumtype) rem) << (8 * i); \ 110 retval |= ((accumtype) rem) << (8 * i); \
98 } \ 111 } \
99 MP_MFREE (quo); \ 112 bignum_fini (quo); \
100 } while (0) 113 } while (0)
101 114
102 int 115 int
103 bignum_to_int (bignum b) 116 bignum_to_int (bignum b)
104 { 117 {
105 short rem, sign; 118 short rem, sign;
106 unsigned int retval = 0; 119 unsigned int retval = 0;
107 REGISTER unsigned int i; 120 REGISTER unsigned int i;
108 MINT *quo; 121 MINT *quo;
109 122
110 sign = MP_MCMP (b, bignum_zero) < 0 ? -1 : 1; 123 sign = bignum_sign (b);
111 quo = MP_ITOM (sign);
112 BIGNUM_TO_TYPE (int, unsigned int); 124 BIGNUM_TO_TYPE (int, unsigned int);
113 return ((int) retval) * sign; 125 return ((int) retval) * sign;
114 } 126 }
115 127
116 unsigned int 128 unsigned int
117 bignum_to_uint (bignum b) 129 bignum_to_uint (bignum b)
118 { 130 {
119 short rem; 131 short rem, sign;
120 unsigned int retval = 0U; 132 unsigned int retval = 0U;
121 REGISTER unsigned int i; 133 REGISTER unsigned int i;
122 MINT *quo; 134 MINT *quo;
123 135
124 quo = MP_ITOM (MP_MCMP (b, bignum_zero) < 0 ? -1 : 1); 136 sign = bignum_sign (b);
125 BIGNUM_TO_TYPE (unsigned int, unsigned int); 137 BIGNUM_TO_TYPE (unsigned int, unsigned int);
126 return retval; 138 return retval;
127 } 139 }
128 140
129 long 141 long
132 short rem, sign; 144 short rem, sign;
133 unsigned long retval = 0L; 145 unsigned long retval = 0L;
134 REGISTER unsigned int i; 146 REGISTER unsigned int i;
135 MINT *quo; 147 MINT *quo;
136 148
137 sign = MP_MCMP (b, bignum_zero) < 0 ? -1 : 1; 149 sign = bignum_sign (b);
138 quo = MP_ITOM (sign);
139 BIGNUM_TO_TYPE (long, unsigned long); 150 BIGNUM_TO_TYPE (long, unsigned long);
140 return ((long) retval) * sign; 151 return ((long) retval) * sign;
141 } 152 }
142 153
143 unsigned long 154 unsigned long
144 bignum_to_ulong (bignum b) 155 bignum_to_ulong (bignum b)
145 { 156 {
146 short rem; 157 short rem, sign;
147 unsigned long retval = 0UL; 158 unsigned long retval = 0UL;
148 REGISTER unsigned int i; 159 REGISTER unsigned int i;
149 MINT *quo; 160 MINT *quo;
150 161
151 quo = MP_ITOM (MP_MCMP (b, bignum_zero) < 0 ? -1 : 1); 162 sign = bignum_sign (b);
152 BIGNUM_TO_TYPE (unsigned long, unsigned long); 163 BIGNUM_TO_TYPE (unsigned long, unsigned long);
153 return retval; 164 return retval;
154 } 165 }
155 166
156 double 167 double
159 short rem, sign; 170 short rem, sign;
160 double retval = 0.0, factor = 1.0; 171 double retval = 0.0, factor = 1.0;
161 REGISTER unsigned int i; 172 REGISTER unsigned int i;
162 MINT *quo; 173 MINT *quo;
163 174
164 sign = MP_MCMP (b, bignum_zero) < 0 ? -1 : 1; 175 sign = bignum_sign (b);
165 quo = MP_ITOM (sign); 176 bignum_init (quo);
166 MP_MULT (b, quo, quo); 177 if (sign < 0)
178 {
179 MP_MSUB (bignum_zero, b, quo);
180 }
181 else
182 {
183 MP_MOVE (b, quo);
184 }
185
167 for (i = 0U; MP_MCMP (quo, bignum_zero) > 0; i++) 186 for (i = 0U; MP_MCMP (quo, bignum_zero) > 0; i++)
168 { 187 {
169 MP_SDIV (quo, 256, quo, &rem); 188 MP_SDIV (quo, 256, quo, &rem);
170 retval += rem * factor; 189 retval += rem * factor;
171 factor *= 256.0; 190 factor *= 256.0;
301 } 320 }
302 321
303 void bignum_ceil (bignum quotient, bignum N, bignum D) 322 void bignum_ceil (bignum quotient, bignum N, bignum D)
304 { 323 {
305 MP_MDIV (N, D, quotient, intern_bignum); 324 MP_MDIV (N, D, quotient, intern_bignum);
306 if (MP_MCMP (intern_bignum, bignum_zero) > 0 && 325 MP_MDIV (N, D, quotient, intern_bignum);
307 MP_MCMP (quotient, bignum_zero) > 0) 326 if (MP_MCMP (intern_bignum, bignum_zero) != 0)
308 MP_MADD (quotient, bignum_one, quotient); 327 {
328 short signN = MP_MCMP (N, bignum_zero);
329 short signD = MP_MCMP (D, bignum_zero);
330
331 /* If the quotient is positive, add one, since we're rounding to
332 positive infinity. */
333 if (signD < 0)
334 {
335 if (signN <= 0)
336 {
337 MP_MADD (quotient, bignum_one, quotient);
338 }
339 }
340 else if (signN >= 0)
341 {
342 MP_MADD (quotient, bignum_one, quotient);
343 }
344 }
309 } 345 }
310 346
311 void bignum_floor (bignum quotient, bignum N, bignum D) 347 void bignum_floor (bignum quotient, bignum N, bignum D)
312 { 348 {
313 MP_MDIV (N, D, quotient, intern_bignum); 349 MP_MDIV (N, D, quotient, intern_bignum);
314 if (MP_MCMP (intern_bignum, bignum_zero) > 0 && 350
315 MP_MCMP (quotient, bignum_zero) < 0) 351 if (MP_MCMP (intern_bignum, bignum_zero) != 0)
316 MP_MSUB (quotient, bignum_one, quotient); 352 {
353 short signN = MP_MCMP (N, bignum_zero);
354 short signD = MP_MCMP (D, bignum_zero);
355
356 /* If the quotient is negative, subtract one, we're rounding to minus
357 infinity. */
358 if (signD < 0)
359 {
360 if (signN >= 0)
361 {
362 MP_MSUB (quotient, bignum_one, quotient);
363 }
364 }
365 else if (signN < 0)
366 {
367 MP_MSUB (quotient, bignum_one, quotient);
368 }
369 }
317 } 370 }
318 371
319 /* RESULT = N to the POWth power */ 372 /* RESULT = N to the POWth power */
320 void 373 void
321 bignum_pow (bignum result, bignum n, unsigned long pow) 374 bignum_pow (bignum result, bignum n, unsigned long pow)