1 // arith11.cpp Copyright (C) 1990-2021 Codemist
2
3 //
4 // Arithmetic functions.
5 // remainder, =,
6 // minusp, plusp
7 //
8 //
9
10 /**************************************************************************
11 * Copyright (C) 2021, Codemist. A C Norman *
12 * *
13 * Redistribution and use in source and binary forms, with or without *
14 * modification, are permitted provided that the following conditions are *
15 * met: *
16 * *
17 * * Redistributions of source code must retain the relevant *
18 * copyright notice, this list of conditions and the following *
19 * disclaimer. *
20 * * Redistributions in binary form must reproduce the above *
21 * copyright notice, this list of conditions and the following *
22 * disclaimer in the documentation and/or other materials provided *
23 * with the distribution. *
24 * *
25 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS *
26 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT *
27 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS *
28 * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE *
29 * COPYRIGHT OWNERS OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, *
30 * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, *
31 * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS *
32 * OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND *
33 * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR *
34 * TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF *
35 * THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH *
36 * DAMAGE. *
37 *************************************************************************/
38
39
40 // $Id: arith11.cpp 5598 2021-01-18 17:27:01Z arthurcnorman $
41
42
43 #include "headers.h"
44
45
46
rembi(LispObject a,LispObject b)47 LispObject rembi(LispObject a, LispObject b)
48 { if (b == fixnum_of_int(0)) return aerror2("bad arg for remainder", a, b);
49 else if (b == fixnum_of_int(1) ||
50 b == fixnum_of_int(-1)) return fixnum_of_int(0);
51 intptr_t n = int_of_fixnum(b);
52 if (signed31_in_ptr(n))
53 { quotbn1(a, n);
54 return fixnum_of_int(nwork);
55 }
56 quotbb(a, make_fake_bignum(n), QUOTBB_REMAINDER_NEEDED);
57 return mv_2;
58 }
59
rembb(LispObject a,LispObject b)60 LispObject rembb(LispObject a, LispObject b)
61 { quotbb(a, b, QUOTBB_REMAINDER_NEEDED);
62 return mv_2;
63 }
64
remis(LispObject a,LispObject b)65 static LispObject remis(LispObject a, LispObject b)
66 { return aerror2("bad arg for remainder", a, b);
67 }
68
remir(LispObject a,LispObject b)69 static LispObject remir(LispObject a, LispObject b)
70 { return aerror2("bad arg for remainder", a, b);
71 }
72
remif(LispObject a,LispObject b)73 static LispObject remif(LispObject a, LispObject b)
74 { return aerror2("bad arg for remainder", a, b);
75 }
76
remsi(LispObject a,LispObject b)77 static LispObject remsi(LispObject a, LispObject b)
78 { return aerror2("bad arg for remainder", a, b);
79 }
80
remsb(LispObject a,LispObject b)81 static LispObject remsb(LispObject a, LispObject b)
82 { return aerror2("bad arg for remainder", a, b);
83 }
84
remsr(LispObject a,LispObject b)85 static LispObject remsr(LispObject a, LispObject b)
86 { return aerror2("bad arg for remainder", a, b);
87 }
88
remsf(LispObject a,LispObject b)89 static LispObject remsf(LispObject a, LispObject b)
90 { return aerror2("bad arg for remainder", a, b);
91 }
92
rembs(LispObject a,LispObject b)93 static LispObject rembs(LispObject a, LispObject b)
94 { return aerror2("bad arg for remainder", a, b);
95 }
96
rembr(LispObject a,LispObject b)97 static LispObject rembr(LispObject a, LispObject b)
98 { return aerror2("bad arg for remainder", a, b);
99 }
100
rembf(LispObject a,LispObject b)101 static LispObject rembf(LispObject a, LispObject b)
102 { return aerror2("bad arg for remainder", a, b);
103 }
104
remri(LispObject a,LispObject b)105 static LispObject remri(LispObject a, LispObject b)
106 { return aerror2("bad arg for remainder", a, b);
107 }
108
remrs(LispObject a,LispObject b)109 static LispObject remrs(LispObject a, LispObject b)
110 { return aerror2("bad arg for remainder", a, b);
111 }
112
remrb(LispObject a,LispObject b)113 static LispObject remrb(LispObject a, LispObject b)
114 { return aerror2("bad arg for remainder", a, b);
115 }
116
remrr(LispObject a,LispObject b)117 static LispObject remrr(LispObject a, LispObject b)
118 { return aerror2("bad arg for remainder", a, b);
119 }
120
remrf(LispObject a,LispObject b)121 static LispObject remrf(LispObject a, LispObject b)
122 { return aerror2("bad arg for remainder", a, b);
123 }
124
remfi(LispObject a,LispObject b)125 static LispObject remfi(LispObject a, LispObject b)
126 { return aerror2("bad arg for remainder", a, b);
127 }
128
remfs(LispObject a,LispObject b)129 static LispObject remfs(LispObject a, LispObject b)
130 { return aerror2("bad arg for remainder", a, b);
131 }
132
remfb(LispObject a,LispObject b)133 static LispObject remfb(LispObject a, LispObject b)
134 { return aerror2("bad arg for remainder", a, b);
135 }
136
remfr(LispObject a,LispObject b)137 static LispObject remfr(LispObject a, LispObject b)
138 { return aerror2("bad arg for remainder", a, b);
139 }
140
remff(LispObject a,LispObject b)141 static LispObject remff(LispObject a, LispObject b)
142 { return aerror2("bad arg for remainder", a, b);
143 }
144
Cremainder(LispObject a,LispObject b)145 LispObject Cremainder(LispObject a, LispObject b)
146 { intptr_t c;
147 switch (static_cast<int>(a) & XTAG_BITS)
148 { case TAG_FIXNUM:
149 switch (static_cast<int>(b) & XTAG_BITS)
150 { case TAG_FIXNUM:
151 //
152 // This is where fixnum % fixnum arithmetic happens - the case I most want to
153 // make efficient.
154 //
155 if (b == fixnum_of_int(0))
156 return aerror2("bad arg for remainder", a, b);
157 // No overflow is possible in a remaindering operation
158 { intptr_t aa = int_of_fixnum(a);
159 intptr_t bb = int_of_fixnum(b);
160 c = aa % bb;
161 //
162 // C does not specify just what happens when % is used with negative
163 // operands (except maybe if the division went exactly), so here I do
164 // some adjusting, assuming that the quotient returned was one of the
165 // integral values surrounding the exact result.
166 //
167 if (aa < 0)
168 { if (c > 0) c -= bb;
169 }
170 else if (c < 0) c += bb;
171 return fixnum_of_int(c);
172 }
173 //
174 // Common Lisp defines a meaning for the remainder function when applied
175 // to floating point values - so there is a whole pile of mess here to
176 // support that. Standard Lisp is only concerned with fixnums and
177 // bignums, but can tolerate the extra generality.
178 //
179 case XTAG_SFLOAT:
180 return remis(a, b);
181 case TAG_NUMBERS:
182 case TAG_NUMBERS+TAG_XBIT:
183 { int32_t hb = type_of_header(numhdr(b));
184 switch (hb)
185 { case TYPE_BIGNUM:
186 //
187 // When I divide a fixnum a by a bignum b the remainder is a except in
188 // the case that a = 0xf8000000 and b = 0x08000000 (for 28-bit fixnums)
189 // or a = -2^59 and b = 2^59 on a 64-bit system, in which case the answer
190 // is zero. This is all because the most negative fixnum has magnitude
191 // such that when you negated it you ended up with a bignum.
192 //
193 if (a == MOST_NEGATIVE_FIXNUM &&
194 ((SIXTY_FOUR_BIT &&
195 bignum_length(b) == CELL+8 &&
196 bignum_digits(b)[0] == 0 &&
197 bignum_digits(b)[1] == (1<<(59-31))) ||
198 (!SIXTY_FOUR_BIT &&
199 bignum_length(b) == CELL+4 &&
200 bignum_digits(b)[0] == (1<<27))))
201 return fixnum_of_int(0);
202 else return a;
203 case TYPE_RATNUM:
204 return remir(a, b);
205 default:
206 return aerror1("Bad arg for remainder", b);
207 }
208 }
209 case TAG_BOXFLOAT:
210 case TAG_BOXFLOAT+TAG_XBIT:
211 return remif(a, b);
212 default:
213 return aerror1("Bad arg for remainder", b);
214 }
215 case XTAG_SFLOAT:
216 switch (static_cast<int>(b) & XTAG_BITS)
217 { case TAG_FIXNUM:
218 return remsi(a, b);
219 case XTAG_SFLOAT:
220 { LispObject q = pack_immediate_float(
221 value_of_immediate_float(a) /
222 value_of_immediate_float(b), a, b);
223 return pack_immediate_float(
224 value_of_immediate_float(a) -
225 (value_of_immediate_float(b) *
226 value_of_immediate_float(q)), a, b);
227 }
228 case TAG_NUMBERS:
229 case TAG_NUMBERS+TAG_XBIT:
230 { int32_t hb = type_of_header(numhdr(b));
231 switch (hb)
232 { case TYPE_BIGNUM:
233 return remsb(a, b);
234 case TYPE_RATNUM:
235 return remsr(a, b);
236 default:
237 return aerror1("Bad arg for remainder", b);
238 }
239 }
240 case TAG_BOXFLOAT:
241 case TAG_BOXFLOAT+TAG_XBIT:
242 return remsf(a, b);
243 default:
244 return aerror1("Bad arg for remainder", b);
245 }
246 case TAG_NUMBERS:
247 case TAG_NUMBERS+TAG_XBIT:
248 { int32_t ha = type_of_header(numhdr(a));
249 switch (ha)
250 { case TYPE_BIGNUM:
251 switch (static_cast<int>(b) & XTAG_BITS)
252 { case TAG_FIXNUM:
253 return rembi(a, b);
254 case XTAG_SFLOAT:
255 return rembs(a, b);
256 case TAG_NUMBERS:
257 case TAG_NUMBERS+TAG_XBIT:
258 { int32_t hb = type_of_header(numhdr(b));
259 switch (hb)
260 { case TYPE_BIGNUM:
261 return rembb(a, b);
262 case TYPE_RATNUM:
263 return rembr(a, b);
264 default:
265 return aerror1("Bad arg for remainder", b);
266 }
267 }
268 case TAG_BOXFLOAT:
269 case TAG_BOXFLOAT+TAG_XBIT:
270 return rembf(a, b);
271 default:
272 return aerror1("Bad arg for remainder", b);
273 }
274 case TYPE_RATNUM:
275 switch (static_cast<int>(b) & XTAG_BITS)
276 { case TAG_FIXNUM:
277 return remri(a, b);
278 case XTAG_SFLOAT:
279 return remrs(a, b);
280 case TAG_NUMBERS:
281 case TAG_NUMBERS+TAG_XBIT:
282 { int32_t hb = type_of_header(numhdr(b));
283 switch (hb)
284 { case TYPE_BIGNUM:
285 return remrb(a, b);
286 case TYPE_RATNUM:
287 return remrr(a, b);
288 default:
289 return aerror1("Bad arg for remainder", b);
290 }
291 }
292 case TAG_BOXFLOAT:
293 case TAG_BOXFLOAT+TAG_XBIT:
294 return remrf(a, b);
295 default:
296 return aerror1("Bad arg for remainder", b);
297 }
298 default: return aerror1("Bad arg for remainder", a);
299 }
300 }
301 case TAG_BOXFLOAT:
302 case TAG_BOXFLOAT+TAG_XBIT:
303 switch (static_cast<int>(b) & TAG_BITS)
304 { case TAG_FIXNUM:
305 return remfi(a, b);
306 case XTAG_SFLOAT:
307 return remfs(a, b);
308 case TAG_NUMBERS:
309 case TAG_NUMBERS+TAG_XBIT:
310 { int32_t hb = type_of_header(numhdr(b));
311 switch (hb)
312 { case TYPE_BIGNUM:
313 return remfb(a, b);
314 case TYPE_RATNUM:
315 return remfr(a, b);
316 default:
317 return aerror1("Bad arg for remainder", b);
318 }
319 }
320 case TAG_BOXFLOAT:
321 case TAG_BOXFLOAT+TAG_XBIT:
322 return remff(a, b);
323 default:
324 return aerror1("Bad arg for remainder", b);
325 }
326 default:
327 return aerror1("Bad arg for remainder", a);
328 }
329 }
330
331 //
332 // In the cases that I expect to be most speed-critical I will
333 // implement "mod" directly. But in a load of other cases I will just
334 // activate the existing "remainder" code and then make a few final
335 // adjustments. This MAY lead to error messages (on modulus by zero)
336 // mentioning remainder rather than mod....
337 // I will leave in the whole structure of separate functions for each
338 // case since that will be useful if I ever need to come back here and
339 // fine-tune more of the type-combinations. As a first pass I give
340 // special treatment to (fixnum,fixnum) and (bignum,fixnum)
341 //
342
mod_by_rem(LispObject a,LispObject b)343 static LispObject mod_by_rem(LispObject a, LispObject b)
344 { bool sb = minusp(b);
345 { Save save(b);
346 a = Cremainder(a, b); // Repeats dispatch on argument type. Sorry
347 errexit();
348 save.restore(b);
349 }
350 if (sb)
351 { if (plusp(a)) a = plus2(a, b);
352 }
353 else if (minusp(a)) a = plus2(a, b);
354 return a;
355 }
356
modib(LispObject a,LispObject b)357 static LispObject modib(LispObject a, LispObject b)
358 { return mod_by_rem(a, b);
359 }
360
modbi(LispObject a,LispObject b)361 static LispObject modbi(LispObject a, LispObject b)
362 { if (b == fixnum_of_int(0)) return aerror2("bad arg for mod", a, b);
363 if (b == fixnum_of_int(1) || b == fixnum_of_int(-1))
364 return fixnum_of_int(0);
365 intptr_t n = int_of_fixnum(b);
366 if (signed31_in_ptr(n))
367 { quotbn1(a, n);
368 // The modulus must have the same sign as b (ie as n).
369 if (n < 0)
370 { if (nwork > 0) nwork += n;
371 }
372 else if (nwork < 0) nwork += n;
373 return fixnum_of_int(nwork);
374 }
375 quotbb(a, make_fake_bignum(n), QUOTBB_REMAINDER_NEEDED);
376 // here b had been a fixnum and the absolute value of the remainder in mv_2
377 // must be strictly smaller than |b|. And if I need to do and adjustment to
378 // turn remainder into modulus that will still leave a value that will be
379 // smaller that b was, and hence it must be a fixnum. That means I can do
380 // the adjustment arithmetic cheaply!
381 if (n < 0)
382 { if (plusp(mv_2)) mv_2 = mv_2 + b - TAG_FIXNUM;
383 }
384 else if (minusp(mv_2)) mv_2 = mv_2 + b - TAG_FIXNUM;
385 return mv_2;
386 }
387
modbb(LispObject a,LispObject b)388 static LispObject modbb(LispObject a, LispObject b)
389 { return mod_by_rem(a, b);
390 }
391
modis(LispObject a,LispObject b)392 static LispObject modis(LispObject a, LispObject b)
393 { return mod_by_rem(a, b);
394 }
395
modir(LispObject a,LispObject b)396 static LispObject modir(LispObject a, LispObject b)
397 { return mod_by_rem(a, b);
398 }
399
modif(LispObject a,LispObject b)400 static LispObject modif(LispObject a, LispObject b)
401 { return mod_by_rem(a, b);
402 }
403
modsi(LispObject a,LispObject b)404 static LispObject modsi(LispObject a, LispObject b)
405 { return mod_by_rem(a, b);
406 }
407
modsb(LispObject a,LispObject b)408 static LispObject modsb(LispObject a, LispObject b)
409 { return mod_by_rem(a, b);
410 }
411
modsr(LispObject a,LispObject b)412 static LispObject modsr(LispObject a, LispObject b)
413 { return mod_by_rem(a, b);
414 }
415
modsf(LispObject a,LispObject b)416 static LispObject modsf(LispObject a, LispObject b)
417 { return mod_by_rem(a, b);
418 }
419
modbs(LispObject a,LispObject b)420 static LispObject modbs(LispObject a, LispObject b)
421 { return mod_by_rem(a, b);
422 }
423
modbr(LispObject a,LispObject b)424 static LispObject modbr(LispObject a, LispObject b)
425 { return mod_by_rem(a, b);
426 }
427
modbf(LispObject a,LispObject b)428 static LispObject modbf(LispObject a, LispObject b)
429 { return mod_by_rem(a, b);
430 }
431
modri(LispObject a,LispObject b)432 static LispObject modri(LispObject a, LispObject b)
433 { return mod_by_rem(a, b);
434 }
435
modrs(LispObject a,LispObject b)436 static LispObject modrs(LispObject a, LispObject b)
437 { return mod_by_rem(a, b);
438 }
439
modrb(LispObject a,LispObject b)440 static LispObject modrb(LispObject a, LispObject b)
441 { return mod_by_rem(a, b);
442 }
443
modrr(LispObject a,LispObject b)444 static LispObject modrr(LispObject a, LispObject b)
445 { return mod_by_rem(a, b);
446 }
447
modrf(LispObject a,LispObject b)448 static LispObject modrf(LispObject a, LispObject b)
449 { return mod_by_rem(a, b);
450 }
451
modfi(LispObject a,LispObject b)452 static LispObject modfi(LispObject a, LispObject b)
453 { return mod_by_rem(a, b);
454 }
455
modfs(LispObject a,LispObject b)456 static LispObject modfs(LispObject a, LispObject b)
457 { return mod_by_rem(a, b);
458 }
459
modfb(LispObject a,LispObject b)460 static LispObject modfb(LispObject a, LispObject b)
461 { return mod_by_rem(a, b);
462 }
463
modfr(LispObject a,LispObject b)464 static LispObject modfr(LispObject a, LispObject b)
465 { return mod_by_rem(a, b);
466 }
467
ccl_modff(LispObject a,LispObject b)468 static LispObject ccl_modff(LispObject a, LispObject b)
469 { return mod_by_rem(a, b);
470 }
471
modulus(LispObject a,LispObject b)472 LispObject modulus(LispObject a, LispObject b)
473 { switch (static_cast<int>(a) & XTAG_BITS)
474 { case TAG_FIXNUM:
475 switch (static_cast<int>(b) & XTAG_BITS)
476 { case TAG_FIXNUM:
477 //
478 // This is where fixnum % fixnum arithmetic happens - the case I most want to
479 // make efficient.
480 //
481 { intptr_t p = int_of_fixnum(a);
482 intptr_t q = int_of_fixnum(b);
483 if (q == 0) return aerror2("bad arg for mod", a, b);
484 p = p % q;
485 if (q < 0)
486 { if (p > 0) p += q;
487 }
488 else if (p < 0) p += q;
489 // No overflow is possible in a modulus operation
490 return fixnum_of_int(p);
491 }
492 //
493 // Common Lisp defines a meaning for the modulus function when applied
494 // to floating point values - so there is a whole pile of mess here to
495 // support that. Standard Lisp is only concerned with fixnums and
496 // bignums.
497 //
498 case XTAG_SFLOAT:
499 return modis(a, b);
500 case TAG_NUMBERS:
501 case TAG_NUMBERS+TAG_XBIT:
502 { int hb = type_of_header(numhdr(b));
503 switch (hb)
504 { case TYPE_BIGNUM:
505 return modib(a, b);
506 case TYPE_RATNUM:
507 return modir(a, b);
508 default:
509 return aerror1("Bad arg for mod", b);
510 }
511 }
512 case TAG_BOXFLOAT:
513 case TAG_BOXFLOAT+TAG_XBIT:
514 return modif(a, b);
515 default:
516 return aerror1("Bad arg for mod", b);
517 }
518 case XTAG_SFLOAT:
519 switch (static_cast<int>(b) & XTAG_BITS)
520 { case TAG_FIXNUM:
521 return modsi(a, b);
522 case XTAG_SFLOAT:
523 { LispObject q = pack_immediate_float(
524 value_of_immediate_float(a) /
525 value_of_immediate_float(b), a, b);
526 // What I have here is a remainder rather than a modulus. I ought to adjust
527 // this code in a spirit of directed rounding.
528 double r = value_of_immediate_float(a) -
529 value_of_immediate_float(b) *
530 value_of_immediate_float(q);
531 return pack_immediate_float(r, a, b);
532 }
533 case TAG_NUMBERS:
534 case TAG_NUMBERS+TAG_XBIT:
535 { int32_t hb = type_of_header(numhdr(b));
536 switch (hb)
537 { case TYPE_BIGNUM:
538 return modsb(a, b);
539 case TYPE_RATNUM:
540 return modsr(a, b);
541 default:
542 return aerror1("Bad arg for mod", b);
543 }
544 }
545 case TAG_BOXFLOAT:
546 case TAG_BOXFLOAT+TAG_XBIT:
547 return modsf(a, b);
548 default:
549 return aerror1("Bad arg for mod", b);
550 }
551 case TAG_NUMBERS:
552 case TAG_NUMBERS+TAG_XBIT:
553 { int ha = type_of_header(numhdr(a));
554 switch (ha)
555 { case TYPE_BIGNUM:
556 switch (static_cast<int>(b) & XTAG_BITS)
557 { case TAG_FIXNUM:
558 return modbi(a, b);
559 case XTAG_SFLOAT:
560 return modbs(a, b);
561 case TAG_NUMBERS:
562 case TAG_NUMBERS+TAG_XBIT:
563 { int hb = type_of_header(numhdr(b));
564 switch (hb)
565 { case TYPE_BIGNUM:
566 return modbb(a, b);
567 case TYPE_RATNUM:
568 return modbr(a, b);
569 default:
570 return aerror1("Bad arg for mod", b);
571 }
572 }
573 case TAG_BOXFLOAT:
574 case TAG_BOXFLOAT+TAG_XBIT:
575 return modbf(a, b);
576 default:
577 return aerror1("Bad arg for mod", b);
578 }
579 case TYPE_RATNUM:
580 switch (static_cast<int>(b) & XTAG_BITS)
581 { case TAG_FIXNUM:
582 return modri(a, b);
583 case XTAG_SFLOAT:
584 return modrs(a, b);
585 case TAG_NUMBERS:
586 case TAG_NUMBERS+TAG_XBIT:
587 { int hb = type_of_header(numhdr(b));
588 switch (hb)
589 { case TYPE_BIGNUM:
590 return modrb(a, b);
591 case TYPE_RATNUM:
592 return modrr(a, b);
593 default:
594 return aerror1("Bad arg for mod", b);
595 }
596 }
597 case TAG_BOXFLOAT:
598 case TAG_BOXFLOAT+TAG_XBIT:
599 return modrf(a, b);
600 default:
601 return aerror1("Bad arg for mod", b);
602 }
603 default: return aerror1("Bad arg for mod", a);
604 }
605 }
606 case TAG_BOXFLOAT:
607 case TAG_BOXFLOAT+TAG_XBIT:
608 switch (static_cast<int>(b) & XTAG_BITS)
609 { case TAG_FIXNUM:
610 return modfi(a, b);
611 case XTAG_SFLOAT:
612 return modfs(a, b);
613 case TAG_NUMBERS:
614 case TAG_NUMBERS+TAG_XBIT:
615 { int hb = type_of_header(numhdr(b));
616 switch (hb)
617 { case TYPE_BIGNUM:
618 return modfb(a, b);
619 case TYPE_RATNUM:
620 return modfr(a, b);
621 default:
622 return aerror1("Bad arg for mod", b);
623 }
624 }
625 case TAG_BOXFLOAT:
626 case TAG_BOXFLOAT+TAG_XBIT:
627 return ccl_modff(a, b);
628 default:
629 return aerror1("Bad arg for mod", b);
630 }
631 default:
632 return aerror1("Bad arg for mod", a);
633 }
634 }
635
zerop(LispObject a)636 bool zerop(LispObject a)
637 { switch (static_cast<int>(a) & XTAG_BITS)
638 { case TAG_FIXNUM:
639 return (a == fixnum_of_int(0));
640 case TAG_NUMBERS:
641 case TAG_NUMBERS+TAG_XBIT:
642 // #C(r i) must satisfy zerop is r and i both do
643 if (is_complex(a) && zerop(real_part(a)))
644 return zerop(imag_part(a));
645 else return false;
646 case XTAG_SFLOAT:
647 return value_of_immediate_float(a) == 0.0;
648 case TAG_BOXFLOAT:
649 case TAG_BOXFLOAT+TAG_XBIT:
650 return float_of_number(a) == 0.0;
651 default:
652 return false;
653 }
654 }
655
onep(LispObject a)656 bool onep(LispObject a)
657 { switch (static_cast<int>(a) & XTAG_BITS)
658 { case TAG_FIXNUM:
659 return (a == fixnum_of_int(1));
660 case TAG_NUMBERS:
661 case TAG_NUMBERS+TAG_XBIT:
662 // #C(r i) must satisfy onep(r) and zerop(i)
663 if (is_complex(a) && onep(real_part(a)))
664 return zerop(imag_part(a));
665 else return false;
666 case XTAG_SFLOAT:
667 return value_of_immediate_float(a) == 1.0;
668 case TAG_BOXFLOAT:
669 case TAG_BOXFLOAT+TAG_XBIT:
670 return (float_of_number(a) == 1.0);
671 default:
672 return false;
673 }
674 }
675
676 //
677 // sign testing
678 //
679
minusp(LispObject a)680 bool minusp(LispObject a)
681 { switch (static_cast<int>(a) & XTAG_BITS)
682 { case TAG_FIXNUM:
683 return fixnum_minusp(a);
684 case XTAG_SFLOAT:
685 return value_of_immediate_float(a) < 0.0;
686 case TAG_NUMBERS:
687 case TAG_NUMBERS+TAG_XBIT:
688 { int ha = type_of_header(numhdr(a));
689 switch (ha)
690 { case TYPE_BIGNUM:
691 return bignum_minusp(a);
692 case TYPE_RATNUM:
693 return minusp(numerator(a));
694 default:
695 return aerror1("Bad arg for minusp", a);
696 return 0;
697 }
698 }
699 case TAG_BOXFLOAT:
700 case TAG_BOXFLOAT+TAG_XBIT:
701 return float_of_number(a) < 0.0;
702 default:
703 return aerror1("Bad arg for minusp", a);
704 return 0;
705 }
706 }
707
708
plusp(LispObject a)709 bool plusp(LispObject a)
710 { switch (static_cast<int>(a) & XTAG_BITS)
711 { case TAG_FIXNUM:
712 return (intptr_t)a > fixnum_of_int(0);
713 case XTAG_SFLOAT:
714 return value_of_immediate_float(a) > 0.0;
715 case TAG_NUMBERS:
716 case TAG_NUMBERS+TAG_XBIT:
717 { int ha = type_of_header(numhdr(a));
718 switch (ha)
719 { case TYPE_BIGNUM:
720 { size_t l = (bignum_length(a)-CELL-4)/4;
721 // This is OK because a bignum can never have the value zero
722 return ((int32_t)bignum_digits(a)[l] >= 0);
723 }
724 case TYPE_RATNUM:
725 return plusp(numerator(a));
726 default:
727 return aerror1("Bad arg for plusp", a);
728 return 0;
729 }
730 }
731 case TAG_BOXFLOAT:
732 case TAG_BOXFLOAT+TAG_XBIT:
733 return float_of_number(a) > 0.0;
734 default:
735 return aerror1("Bad arg for plusp", a);
736 return 0;
737 }
738 }
739
740
741 //
742 // Numeric equality - note that comparisons involving non-numbers
743 // are errors here (unlike the position in eql, equal, equalp). Also
744 // this must be coded so that it never provokes garbage collection.
745 //
746
numeqis(LispObject a,LispObject b)747 static bool numeqis(LispObject a, LispObject b)
748 { return eq_i64d(int_of_fixnum(a), value_of_immediate_float(b));
749 }
750
numeqic(LispObject a,LispObject b)751 static bool numeqic(LispObject a, LispObject b)
752 { if (!zerop(imag_part(b))) return false;
753 else return numeq2(a, real_part(b));
754 }
755
756 #define numeqif(a,b) eq_i64d(int_of_fixnum(a), float_of_number(b))
757
758 #define numeqsi(a, b) numeqis(b, a)
759
numeqsb(LispObject a,LispObject b)760 static bool numeqsb(LispObject a, LispObject b)
761 //
762 // This is coded to allow comparison of any floating type
763 // with a bignum. Well when I say "any floating type" I think that the
764 // 128 bit one will need further work!
765 //
766 { double d = float_of_number(a), d1;
767 int x;
768 int32_t w;
769 size_t len;
770 uint32_t u;
771 // MOST_NEGATIVE_FIXVAL will be a power of 2, and so even if is bigger
772 // then 2^53 it can convert exactly to a double, or indeed to any floating
773 // point type. My reasoning here is that b is a BIGNUM and hence its
774 // value is outside the range of fixnums. So if the value of a is within
775 // the range of fixnums we can not have equality.
776 if (d >= static_cast<double>(MOST_NEGATIVE_FIXVAL) &&
777 d < -static_cast<double>(MOST_NEGATIVE_FIXVAL)) return false;
778 len = (bignum_length(b)-CELL-4)/4;
779 if (len == 0) // One word bignums can be treated specially
780 { int32_t v = bignum_digits(b)[0];
781 return (d == static_cast<double>(v));
782 }
783 d1 = std::frexp(d, &x); // separate exponent from mantissa
784 if (d1 == 1.0) d1 = 0.5, x++; // For Zortech
785 // The exponent x must be positive here, hence the % operation is defined
786 d1 = std::ldexp(d1, x % 31);
787 //
788 // At most 3 words in the bignum may contain nonzero data - I subtract
789 // the (double) value of those bits off and check that (a) the floating
790 // result left is zero and (b) there are no more bits left.
791 //
792 x = x / 31;
793 if (x < 0 || (size_t)x != len) return false;
794 w = bignum_digits(b)[len];
795 d1 = (d1 - static_cast<double>(w)) * TWO_31;
796 u = bignum_digits(b)[--len];
797 d1 = (d1 - static_cast<double>(u)) * TWO_31;
798 if (len > 0)
799 { u = bignum_digits(b)[--len];
800 d1 = d1 - static_cast<double>(u);
801 }
802 if (d1 != 0.0) return false;
803 while (len != 0)
804 if (bignum_digits(b)[--len] != 0) return false;
805 return true;
806 }
807
numeqsr(LispObject a,LispObject b)808 static bool numeqsr(LispObject a, LispObject b)
809 //
810 // Here I will rely somewhat on the use of IEEE floating point values
811 // (an in particular the weaker supposition that I have floating point
812 // with a binary radix). Then for equality the denominator of b must
813 // be a power of 2, which I can test for and then account for.
814 //
815 { LispObject nb = numerator(b), db = denominator(b);
816 double d = float_of_number(a), d1;
817 int x;
818 int32_t dx, w, len;
819 uint32_t u, bit;
820 //
821 // first I will check that db (which will be positive) is a power of 2,
822 // and set dx to indicate what power of two it is.
823 // Note that db != 0 and that one of the top two words of a bignum
824 // must be nonzero (for normalisation) so I end up with a nonzero
825 // value in the variable 'bit'
826 //
827 if (is_fixnum(db))
828 { bit = int_of_fixnum(db);
829 w = bit;
830 if (w != (w & (-w))) return false; // not a power of 2
831 dx = 0;
832 }
833 else if (is_numbers(db) && is_bignum(db))
834 { int32_t lenb = (bignum_length(db)-CELL-4)/4;
835 bit = bignum_digits(db)[lenb];
836 //
837 // I need to cope with bignums where the leading digits is zero because
838 // the 0x80000000 bit of the next word down is 1. To do this I treat
839 // the number as having one fewer digits.
840 //
841 if (bit == 0) bit = bignum_digits(db)[--lenb];
842 w = bit;
843 if (w != (w & (-w))) return false; // not a power of 2
844 dx = 31*lenb;
845 while (--lenb != 0) // check that the rest of db is zero
846 if (bignum_digits(db)[--lenb] != 0) return false;
847 }
848 else return false; // Odd - what type IS db here? Maybe error.
849 if ((bit & 0xffffU) == 0) dx += 16, bit = bit >> 16;
850 if ((bit & 0xff) == 0) dx += 8, bit = bit >> 8;
851 if ((bit & 0xf) == 0) dx += 4, bit = bit >> 4;
852 if ((bit & 0x3) == 0) dx += 2, bit = bit >> 2;
853 if ((bit & 0x1) == 0) dx += 1;
854 if (is_fixnum(nb))
855 { double d1 = static_cast<double>(int_of_fixnum(nb));
856 //
857 // The ldexp on the next line could potentially underflow. In that case C
858 // defines that the result 0.0 be returned. To avoid trouble I put in a
859 // special test the relies on that fact that a value represented as a rational
860 // would not have been zero.
861 //
862 if (dx > 10000) return false; // Avoid gross underflow
863 d1 = std::ldexp(d1, static_cast<int>(-dx));
864 return (d == d1 && d != 0.0);
865 }
866 len = (bignum_length(nb)-CELL-4)/4;
867 if (len == 0) // One word bignums can be treated specially
868 { int32_t v = bignum_digits(nb)[0];
869 double d1;
870 if (dx > 10000) return false; // Avoid gross underflow
871 d1 = std::ldexp(static_cast<double>(v), static_cast<int>(-dx));
872 return (d == d1 && d != 0.0);
873 }
874 d1 = std::frexp(d, &x); // separate exponent from mantissa
875 if (d1 == 1.0) d1 = 0.5, x++; // For Zortech
876 dx += x; // adjust to allow for the denominator
877 d1 = std::ldexp(d1, static_cast<int>(dx % 31));
878 // can neither underflow nor overflow here
879 //
880 // At most 3 words in the bignum may contain nonzero data - I subtract
881 // the (double) value of those bits off and check that (a) the floating
882 // result left is zero and (b) there are no more bits left.
883 //
884 dx = dx / 31;
885 if (dx != len) return false;
886 w = bignum_digits(nb)[len];
887 d1 = (d1 - static_cast<double>(w)) * TWO_31;
888 u = bignum_digits(nb)[--len];
889 d1 = (d1 - static_cast<double>(u)) * TWO_31;
890 if (len > 0)
891 { u = bignum_digits(nb)[--len];
892 d1 = d1 - static_cast<double>(u);
893 }
894 if (d1 != 0.0) return false;
895 while (len != 0)
896 if (bignum_digits(nb)[--len] != 0) return false;
897 return bignum_digits(nb)[0] == 0;
898 }
899
900 #define numeqsc(a, b) numeqic(a, b)
901
numeqsf(LispObject a,LispObject b)902 static bool numeqsf(LispObject a, LispObject b)
903 { return value_of_immediate_float(a) == float_of_number(b);
904 }
905
906 #define numeqbs(a, b) numeqsb(b, a)
907
numeqbb(LispObject a,LispObject b)908 static bool numeqbb(LispObject a, LispObject b)
909 { int32_t la = bignum_length(a);
910 if (la != (int32_t)bignum_length(b)) return false;
911 la = (la-CELL-4)/4;
912 while (la != 0)
913 if (bignum_digits(a)[la] != bignum_digits(b)[la]) return false;
914 else la--;
915 return bignum_digits(a)[0] == bignum_digits(b)[0];
916 }
917
918 #define numeqbc(a, b) numeqic(a, b)
919
920 #define numeqbf(a, b) numeqsb(b, a)
921
922 #define numeqrs(a, b) numeqsr(b, a)
923
numeqrr(LispObject a,LispObject b)924 static bool numeqrr(LispObject a, LispObject b)
925 { return numeq2(numerator(a), numerator(b)) &&
926 numeq2(denominator(a), denominator(b));
927 }
928
929 #define numeqrc(a, b) numeqic(a, b)
930
931 #define numeqrf(a, b) numeqsr(b, a)
932
933 #define numeqci(a, b) numeqic(b, a)
934
935 #define numeqcs(a, b) numeqic(b, a)
936
937 #define numeqcb(a, b) numeqic(b, a)
938
939 #define numeqcr(a, b) numeqic(b, a)
940
numeqcc(LispObject a,LispObject b)941 static bool numeqcc(LispObject a, LispObject b)
942 { return numeq2(real_part(a), real_part(b)) &&
943 numeq2(imag_part(a), imag_part(b));
944 }
945
946 #define numeqcf(a, b) numeqic(b, a)
947
948 #define numeqfi(a, b) numeqif(b, a)
949
950 #define numeqfs(a, b) numeqsf(b, a)
951
952 #define numeqfb(a, b) numeqbf(b, a)
953
954 #define numeqfr(a, b) numeqrf(b, a)
955
956 #define numeqfc(a, b) numeqic(a, b)
957
numeqff(LispObject a,LispObject b)958 static bool numeqff(LispObject a, LispObject b)
959 { return (float_of_number(a) == float_of_number(b));
960 }
961
962 //
963 // This comparison must signal an error on non-numeric operands in
964 // Common Lisp mode, but behave as EQ in CSL mode.
965 //
966
967 #ifdef COMMON
968 # define differenta return aerror1("Bad arg for =", a); return false
969 # define differentb return aerror1("Bad arg for =", b); return false
970 #else
971 # define differenta return false
972 # define differentb return false
973 #endif
974
975 // This is for the Common Lisp (= u v) case and so it is expected to
976 // accept equality across numeric types... but it entitled to raise an
977 // error if an argument is non-numeric. I might be mor egenerous on the
978 // non-numeric cases.
979
numeq2(LispObject a,LispObject b)980 bool numeq2(LispObject a, LispObject b)
981 { switch (static_cast<int>(a) & XTAG_BITS)
982 { case TAG_FIXNUM:
983 switch (static_cast<int>(b) & XTAG_BITS)
984 { case TAG_FIXNUM:
985 return (a == b);
986 case XTAG_SFLOAT:
987 return numeqis(a, b);
988 case TAG_NUMBERS:
989 case TAG_NUMBERS+TAG_XBIT:
990 { int32_t hb = type_of_header(numhdr(b));
991 switch (hb)
992 { case TYPE_BIGNUM:
993 return false; // fixnum can not be equal to a bignum
994 case TYPE_RATNUM:
995 return false;
996 case TYPE_COMPLEX_NUM:
997 return numeqic(a, b); // (= 2 #C(2.0 0.0))? Yuk
998 default:
999 differentb;
1000 }
1001 }
1002 case TAG_BOXFLOAT:
1003 case TAG_BOXFLOAT+TAG_XBIT:
1004 return numeqif(a, b);
1005 default:
1006 differentb;
1007 }
1008 case XTAG_SFLOAT:
1009 switch (static_cast<int>(b) & XTAG_BITS)
1010 { case TAG_FIXNUM:
1011 return numeqsi(a, b);
1012 case XTAG_SFLOAT:
1013 return value_of_immediate_float(a) ==
1014 value_of_immediate_float(b);
1015 case TAG_NUMBERS:
1016 case TAG_NUMBERS+TAG_XBIT:
1017 { int32_t hb = type_of_header(numhdr(b));
1018 switch (hb)
1019 { case TYPE_BIGNUM:
1020 return numeqsb(a, b);
1021 case TYPE_RATNUM:
1022 return numeqsr(a, b);
1023 case TYPE_COMPLEX_NUM:
1024 return numeqsc(a, b);
1025 default:
1026 differentb;
1027 }
1028 }
1029 case TAG_BOXFLOAT:
1030 case TAG_BOXFLOAT+TAG_XBIT:
1031 return numeqsf(a, b);
1032 default:
1033 differentb;
1034 }
1035 case TAG_NUMBERS:
1036 case TAG_NUMBERS+TAG_XBIT:
1037 { int32_t ha = type_of_header(numhdr(a));
1038 switch (ha)
1039 { case TYPE_BIGNUM:
1040 switch (static_cast<int>(b) & XTAG_BITS)
1041 { case TAG_FIXNUM:
1042 return 0;
1043 case XTAG_SFLOAT:
1044 return numeqbs(a, b);
1045 case TAG_NUMBERS:
1046 case TAG_NUMBERS+TAG_XBIT:
1047 { int32_t hb = type_of_header(numhdr(b));
1048 switch (hb)
1049 { case TYPE_BIGNUM:
1050 return numeqbb(a, b);
1051 case TYPE_RATNUM:
1052 return 0;
1053 case TYPE_COMPLEX_NUM:
1054 return numeqbc(a, b);
1055 default:
1056 differentb;
1057 }
1058 }
1059 case TAG_BOXFLOAT:
1060 case TAG_BOXFLOAT+TAG_XBIT:
1061 return numeqbf(a, b);
1062 default:
1063 differentb;
1064 }
1065 case TYPE_RATNUM:
1066 switch (static_cast<int>(b) & XTAG_BITS)
1067 { case TAG_FIXNUM:
1068 return 0;
1069 case XTAG_SFLOAT:
1070 return numeqrs(a, b);
1071 case TAG_NUMBERS:
1072 case TAG_NUMBERS+TAG_XBIT:
1073 { int32_t hb = type_of_header(numhdr(b));
1074 switch (hb)
1075 { case TYPE_BIGNUM:
1076 return 0;
1077 case TYPE_RATNUM:
1078 return numeqrr(a, b);
1079 case TYPE_COMPLEX_NUM:
1080 return numeqrc(a, b);
1081 default:
1082 differentb;
1083 }
1084 }
1085 case TAG_BOXFLOAT:
1086 case TAG_BOXFLOAT+TAG_XBIT:
1087 return numeqrf(a, b);
1088 default:
1089 differentb;
1090 }
1091 case TYPE_COMPLEX_NUM:
1092 switch (static_cast<int>(b) & XTAG_BITS)
1093 { case TAG_FIXNUM:
1094 return numeqci(a, b);
1095 case XTAG_SFLOAT:
1096 return numeqcs(a, b);
1097 case TAG_NUMBERS:
1098 case TAG_NUMBERS+TAG_XBIT:
1099 { int32_t hb = type_of_header(numhdr(b));
1100 switch (hb)
1101 { case TYPE_BIGNUM:
1102 return numeqcb(a, b);
1103 case TYPE_RATNUM:
1104 return numeqcr(a, b);
1105 case TYPE_COMPLEX_NUM:
1106 return numeqcc(a, b);
1107 default:
1108 differentb;
1109 }
1110 }
1111 case TAG_BOXFLOAT:
1112 case TAG_BOXFLOAT+TAG_XBIT:
1113 return numeqcf(a, b);
1114 default:
1115 differentb;
1116 }
1117 default: differenta;
1118 }
1119 }
1120 case TAG_BOXFLOAT:
1121 case TAG_BOXFLOAT+TAG_XBIT:
1122 switch (static_cast<int>(b) & XTAG_BITS)
1123 { case TAG_FIXNUM:
1124 return numeqfi(a, b);
1125 case XTAG_SFLOAT:
1126 return numeqfs(a, b);
1127 case TAG_NUMBERS:
1128 case TAG_NUMBERS+TAG_XBIT:
1129 { int32_t hb = type_of_header(numhdr(b));
1130 switch (hb)
1131 { case TYPE_BIGNUM:
1132 return numeqfb(a, b);
1133 case TYPE_RATNUM:
1134 return numeqfr(a, b);
1135 case TYPE_COMPLEX_NUM:
1136 return numeqfc(a, b);
1137 default:
1138 differentb;
1139 }
1140 }
1141 case TAG_BOXFLOAT:
1142 case TAG_BOXFLOAT+TAG_XBIT:
1143 return numeqff(a, b);
1144 default:
1145 differentb;
1146 }
1147 default:
1148 differenta;
1149 }
1150 }
1151
1152 // numeq2 as above is willing to compare a float with a bignum and try to
1153 // say if they represent the same value. That is like the Common Lisp
1154 // function (= N1 N2). But Standard Lisp's EQN only accepts numbers as equal
1155 // if they have the same type, and that reduces the complication significantly.
1156
SL_numeq2(LispObject a,LispObject b)1157 bool SL_numeq2(LispObject a, LispObject b)
1158 { unsigned int ha;
1159 switch (static_cast<int>(a) & XTAG_BITS)
1160 { case TAG_FIXNUM:
1161 return (a == b);
1162 case XTAG_SFLOAT:
1163 // OK so Standard Lisp does not have a concept of multiple widths of
1164 // floats, but I am going to make EQN demand that values are the same width
1165 // as well as having the same value for them to be EQN. This is consistent with
1166 // them "having the same type".
1167 if ((static_cast<int>(b) & XTAG_BITS) != XTAG_SFLOAT) return false;
1168 return (value_of_immediate_float(a)==value_of_immediate_float(b));
1169 case TAG_NUMBERS:
1170 case TAG_NUMBERS+TAG_XBIT:
1171 if (a == b) return true;
1172 if ((static_cast<int>(b) & TAG_BITS) != TAG_NUMBERS) return false;
1173 ha = type_of_header(numhdr(a));
1174 if (ha != type_of_header(numhdr(b))) return false;
1175 switch (ha)
1176 { case TYPE_BIGNUM:
1177 return numeqbb(a, b);
1178 case TYPE_RATNUM:
1179 return numeqrr(a, b);
1180 case TYPE_COMPLEX_NUM:
1181 return numeqcc(a, b);
1182 }
1183 differenta;
1184 case TAG_BOXFLOAT:
1185 case TAG_BOXFLOAT+TAG_XBIT:
1186 if ((static_cast<int>(b) & TAG_BITS) != TAG_BOXFLOAT) return false;
1187 ha = type_of_header(flthdr(a));
1188 if (ha != type_of_header(flthdr(b))) return false;
1189 return numeqff(a, b);
1190 default:
1191 return (a == b);
1192 }
1193 }
1194
1195 // end of arith11.cpp
1196