1 /****************************************************************/
2 /* file aritaux.c
3
4 ARIBAS interpreter for Arithmetic
5 Copyright (C) 1996 O.Forster
6
7 This program is free software; you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation; either version 2 of the License, or
10 (at your option) any later version.
11
12 This program is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with this program; if not, write to the Free Software
19 Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
20
21 Address of the author
22
23 Otto Forster
24 Math. Institut der LMU
25 Theresienstr. 39
26 D-80333 Muenchen, Germany
27
28 Email forster@rz.mathematik.uni-muenchen.de
29 */
30 /****************************************************************/
31 /*
32 ** aritaux.c
33 ** auxiliary procedures for arithmetic
34 **
35 ** date of last change
36 ** 1994-12-31
37 ** 2000-12-30: multiprec floats
38 ** 2002-02-16: changesign
39 ** 2002-04-21: chkintnz
40 */
41
42 #include "common.h"
43
44 PUBLIC int setfltprec _((int prec));
45 PUBLIC int deffltprec _((void));
46 PUBLIC int maxfltprec _((void));
47 PUBLIC int fltpreccode _((int prec));
48 PUBLIC int fltprec _((int type));
49 PUBLIC int refnumtrunc _((int prec, truc *ptr, numdata *nptr));
50 PUBLIC int getnumtrunc _((int prec, truc *ptr, numdata *nptr));
51 PUBLIC int getnumalign _((int prec, truc *ptr, numdata *nptr));
52 PUBLIC int alignfloat _((int prec, numdata *nptr));
53 PUBLIC int alignfix _((int prec, numdata *nptr));
54 PUBLIC void adjustoffs _((numdata *npt1, numdata *npt2));
55 PUBLIC int normfloat _((int prec, numdata *nptr));
56 PUBLIC int multtrunc _((int prec, numdata *npt1, numdata *npt2,
57 word2 *hilf));
58 PUBLIC int divtrunc _((int prec, numdata *npt1, numdata *npt2,
59 word2 *hilf));
60 PUBLIC int pwrtrunc _((int prec, unsigned base, unsigned a,
61 numdata *nptr, word2 *hilf));
62 PUBLIC int float2bcd _((int places, truc *p, numdata *nptr,
63 word2 *hilf));
64 PUBLIC int roundbcd _((int prec, numdata *nptr));
65 PUBLIC int flodec2bin _((int prec, numdata *nptr, word2 *hilf));
66 PUBLIC void int2numdat _((int x, numdata *nptr));
67 PUBLIC void cpynumdat _((numdata *npt1, numdata *npt2));
68 PUBLIC int numposneg _((truc *ptr));
69 PUBLIC truc wipesign _((truc *ptr));
70 PUBLIC truc changesign _((truc *ptr));
71 PUBLIC long intretr _((truc *ptr));
72 PUBLIC int bigref _((truc *ptr, word2 **xp, int *sp));
73 PUBLIC int bigretr _((truc *ptr, word2 *x, int *sp));
74 PUBLIC int twocretr _((truc *ptr, word2 *x));
75 PUBLIC int and2arr _((word2 *x, int n, word2 *y, int m));
76 PUBLIC int or2arr _((word2 *x, int n, word2 *y, int m));
77 PUBLIC int xor2arr _((word2 *x, int n, word2 *y, int m));
78 PUBLIC int xorbitvec _((word2 *x, int n, word2 *y, int m));
79 PUBLIC long bit_length _((word2 *x, int n));
80 PUBLIC int chkintnz _((truc sym, truc *ptr));
81 PUBLIC int chkints _((truc sym, truc *argptr, int n));
82 PUBLIC int chkint _((truc sym, truc *ptr));
83 PUBLIC int chkintt _((truc sym, truc *ptr));
84 PUBLIC int chknums _((truc sym, truc *argptr, int n));
85 PUBLIC int chknum _((truc sym, truc *ptr));
86 PUBLIC int chkintvec _((truc sym, truc *vptr));
87 PUBLIC int chknumvec _((truc sym, truc *vptr));
88
89 PRIVATE long decdigs _((numdata *nptr));
90 PRIVATE int malzehnhoch _((int prec, numdata *nptr, long d, word2 *hilf));
91 PRIVATE int twocadjust _((word2 *x, int n));
92
93 PUBLIC int FltPrec[20] =
94 {2, 4, 8, 12, 16, 24, 32, 40, 48, 56, 64, 80, 96, 112, 128,
95 160, 192, 224, 256, 320};
96
97 #ifdef FPREC_HIGH
98 #define MAXFLTLEVEL 19
99 #else
100 #define MAXFLTLEVEL 10
101 #endif
102
103 PUBLIC int MaxFltLevel = MAXFLTLEVEL;
104
105 PRIVATE int floatprec = 2;
106
107 /*------------------------------------------------------------*/
setfltprec(prec)108 PUBLIC int setfltprec(prec)
109 int prec;
110 {
111 int k;
112
113 k = 0;
114 while((k < MaxFltLevel) && (prec > FltPrec[k]))
115 k++;
116 floatprec = FltPrec[k];
117
118 return(floatprec);
119 }
120 /*-------------------------------------------------------------*/
fltpreccode(prec)121 PUBLIC int fltpreccode(prec)
122 int prec;
123 {
124 int k;
125
126 k = 0;
127 while((k < MaxFltLevel) && (prec > FltPrec[k]))
128 k++;
129 return(k);
130 }
131 /*-------------------------------------------------------------*/
maxfltprec()132 PUBLIC int maxfltprec()
133 {
134 return(FltPrec[MaxFltLevel]);
135 }
136 /*-------------------------------------------------------------*/
deffltprec()137 PUBLIC int deffltprec()
138 {
139 return(floatprec);
140 }
141 /*------------------------------------------------------------*/
fltprec(type)142 PUBLIC int fltprec(type)
143 int type;
144 {
145 int k;
146
147 if(type >= fFLTOBJ) {
148 k = (type & PRECMASK) >> 1;
149 return FltPrec[k];
150 }
151 else
152 return(floatprec);
153 }
154 /*------------------------------------------------------------*/
155 /*
156 ** Erzeugt Referenz auf Zahl in *p mit prec 16-bit-Stellen (ohne Kopie)
157 ** Fuer die Zahl 0 wird nptr->len = 0 und nptr->expo = MOSTNEGEX
158 */
refnumtrunc(prec,p,nptr)159 PUBLIC int refnumtrunc(prec,p,nptr)
160 int prec; /* should be >= 2 */
161 truc *p;
162 numdata *nptr;
163 {
164 static word2 ddd[2];
165
166 struct bigcell *big;
167 struct floatcell *fl;
168 long ex;
169 int len, diff, pcode;
170 int flg = *FLAGPTR(p);
171
172 if(flg == fFIXNUM) {
173 if(!*WORD2PTR(p))
174 goto zeroexit;
175 nptr->digits = WORD2PTR(p);
176 nptr->sign = *SIGNPTR(p);
177 nptr->len = 1;
178 nptr->expo = 0;
179 }
180 else if(flg == fBIGNUM) {
181 big = (struct bigcell *)TAddress(p);
182 nptr->sign = big->signum;
183 len = big->len;
184 ex = diff = (len > prec ? len - prec : 0);
185 nptr->digits = &(big->digi0) + diff;
186 nptr->len = len - diff;
187 nptr->expo = (ex << 4);
188 }
189 else if(flg >= fFLTOBJ) {
190 if(flg & FLTZEROBIT)
191 goto zeroexit;
192 pcode = (flg & PRECMASK);
193 len = FltPrec[pcode>>1];
194 fl = (struct floatcell *)TAddress(p);
195 nptr->sign = (fl->signum & FSIGNBIT ? MINUSBYTE : 0);
196 diff = (len > prec ? len - prec : 0);
197 nptr->digits = &(fl->digi0) + diff;
198 nptr->len = len - diff;
199 ex = fl->expo;
200 if(flg >= fHUGEFLOAT) {
201 ex <<= 7;
202 ex += (fl->signum & HUGEMASK);
203 }
204 nptr->expo = ex + (diff << 4);
205 }
206
207 else {
208 error(voidsym,err_case,voidsym);
209 return(0);
210 }
211 return(nptr->len);
212 zeroexit:
213 nptr->digits = ddd;
214 nptr->len = 0;
215 nptr->sign = 0;
216 nptr->expo = MOSTNEGEX;
217 return(0);
218 }
219 /*------------------------------------------------------------------*/
220 /*
221 ** holt Zahl aus *ptr nach *nptr mit prec 16-bit-Stellen (mit Kopie)
222 */
getnumtrunc(prec,ptr,nptr)223 PUBLIC int getnumtrunc(prec,ptr,nptr)
224 int prec;
225 truc *ptr;
226 numdata *nptr;
227 {
228 int len;
229 word2 *saveptr;
230
231 saveptr = nptr->digits;
232 len = refnumtrunc(prec,ptr,nptr);
233 cpyarr(nptr->digits,len,saveptr);
234 nptr->digits = saveptr;
235 return(len);
236 }
237 /*------------------------------------------------------------------*/
238 /*
239 ** holt Zahl aus *p und normalisiert sie
240 */
getnumalign(prec,p,nptr)241 PUBLIC int getnumalign(prec,p,nptr)
242 int prec;
243 truc *p;
244 numdata *nptr;
245 {
246 getnumtrunc(prec,p,nptr);
247 return(alignfloat(prec,nptr));
248 }
249 /*-------------------------------------------------------------------*/
250 /*
251 ** Die in nptr gegebene Zahl wird so normalisiert, dass nptr->expo
252 ** durch 16 teilbar wird.
253 ** Es werden genau prec 16-bit-Stellen benuetzt, falls Zahl /= 0
254 ** Fuer die Zahl 0 ist der Rueckgabewert 0
255 */
alignfloat(prec,nptr)256 PUBLIC int alignfloat(prec,nptr)
257 int prec;
258 numdata *nptr;
259 {
260 int len, b;
261 long sh;
262
263 len = nptr->len;
264 if(len == 0)
265 return(0);
266 b = nptr->expo & 0xF; /* b != 0 only for floats */
267 sh = prec - len;
268 sh <<= 4;
269 sh -= (b ? 16-b : 0);
270 nptr->expo -= sh;
271 nptr->len = lshiftarr(nptr->digits,len,sh);
272 return(nptr->len); /* nptr->len = prec */
273 }
274 /*-----------------------------------------------------------------*/
275 /*
276 ** Die Zahl in nptr wird so normalisiert, dass
277 ** nptr->expo gleich -(prec<<4) wird.
278 ** Falls aber nptr->len groesser als 2*prec wuerde,
279 ** wird nptr unveraendert gelassen und aERROR zurueckgegeben.
280 ** Sonst ist der Rueckgabewert nptr->len
281 */
alignfix(prec,nptr)282 PUBLIC int alignfix(prec,nptr)
283 int prec;
284 numdata *nptr;
285 {
286 int len;
287 long sh;
288
289 if(!(len = nptr->len))
290 return(0);
291 sh = nptr->expo;
292 if(len + (sh>>4) > prec)
293 return(aERROR);
294 sh += prec<<4;
295 nptr->expo -= sh;
296 nptr->len = lshiftarr(nptr->digits,len,sh);
297 return(nptr->len);
298 }
299 /*-------------------------------------------------------------------*/
300 /*
301 ** Die in npt1 und npt2 gegebenen Zahlen werden so normalisiert,
302 ** dass beide expos gleich werden
303 ** Es erfolgt ein Rechtsshift auf npt1->digits oder npt2->digits
304 ** Nur npt1->expo wird auf den korrekten Wert gesetzt
305 ** Es wird vorausgesetzt, dass beide expos durch 16 teilbar sind
306 */
adjustoffs(npt1,npt2)307 PUBLIC void adjustoffs(npt1,npt2)
308 numdata *npt1, *npt2;
309 {
310 int diff;
311
312 if(npt1->len == 0)
313 npt1->expo = npt2->expo;
314 else if(npt2->len == 0)
315 ;
316 else if((diff = (npt1->expo - npt2->expo) >> 4) > 0) {
317 if(diff >= npt2->len)
318 npt2->len = 0;
319 else {
320 npt2->len -= diff;
321 cpyarr(npt2->digits+diff,npt2->len,npt2->digits);
322 }
323 }
324 else if(diff < 0) {
325 npt1->expo = npt2->expo;
326 if(-diff >= npt1->len)
327 npt1->len = 0;
328 else {
329 npt1->len += diff;
330 cpyarr(npt1->digits-diff,npt1->len,npt1->digits);
331 }
332 }
333 }
334 /*------------------------------------------------------------------*/
335 /*
336 ** Die durch nptr gegebene float-Zahl wird normalisiert, so dass
337 ** die mantisse genau (prec*16) Bits lang wird.
338 ** nptr->expo wird der Exponent bzgl. Basis 2
339 ** Rueckgabewert: 0 fuer die Zahl 0, sonst prec
340 ** !!! arbeitet destruktiv auf nptr !!!
341 */
normfloat(prec,nptr)342 PUBLIC int normfloat(prec,nptr)
343 int prec;
344 numdata *nptr;
345 {
346 int len;
347 long sh;
348
349 len = nptr->len;
350 if(len == 0) {
351 nptr->sign = 0;
352 nptr->expo = MOSTNEGEX;
353 return(0);
354 }
355 nptr->len = prec;
356 sh = prec - len + 1;
357 sh <<= 4;
358 sh -= bitlen(nptr->digits[len-1]);
359 lshiftarr(nptr->digits,len,sh);
360 nptr->expo -= sh;
361 return(prec);
362 }
363 /*-------------------------------------------------------------------*/
364 /*
365 ** Die durch npt1 gegebene Zahl wird mit der durch npt2 gegebenen
366 ** Zahl multipliziert. Produkt in npt1
367 ** Die Laenge des Produkts wird auf <= prec 16-bit-Stellen gekuerzt
368 ** Platz hilf muss 3*prec + 4 lang sein
369 */
multtrunc(prec,npt1,npt2,hilf)370 PUBLIC int multtrunc(prec,npt1,npt2,hilf)
371 int prec;
372 numdata *npt1, *npt2;
373 word2 *hilf;
374 {
375 int n, n1, n2;
376 word2 *x, *y, *z;
377 long exponent, diff;
378
379 n1 = npt1->len;
380 n2 = npt2->len;
381 x = npt1->digits;
382 y = npt2->digits;
383
384 diff = n1 - prec;
385 if(diff > 0) {
386 x += diff;
387 n1 -= diff;
388 npt1->expo += diff << 4;
389 }
390 diff = n2 - prec;
391 if(diff > 0) {
392 y += diff;
393 n2 -= diff;
394 npt2->expo += diff << 4;
395 }
396
397 z = hilf + prec + 2;
398 n = multbig(x,n1,y,n2,z,hilf);
399 diff = (n > prec ? n - prec : 0);
400 n -= diff;
401 cpyarr(z+diff,n,npt1->digits);
402 exponent = npt1->expo + npt2->expo + (diff << 4);
403 if(exponent >= maxfltex)
404 return(aERROR);
405 else if(exponent <= -maxfltex) {
406 exponent = MOSTNEGEX;
407 n = 0;
408 }
409 npt1->len = n;
410 npt1->expo = exponent;
411 if(n == 0)
412 npt1->sign = 0;
413 else if((npt1->sign == npt2->sign) || (npt1->sign && npt2->sign))
414 npt1->sign = 0;
415 else
416 npt1->sign = MINUSBYTE;
417 return(n);
418 }
419 /*-------------------------------------------------------------------*/
divtrunc(prec,npt1,npt2,hilf)420 PUBLIC int divtrunc(prec,npt1,npt2,hilf)
421 int prec;
422 numdata *npt1, *npt2;
423 word2 *hilf;
424 {
425 int n, n1, n2;
426 int rlen;
427 long exponent, diff;
428 word2 *x, *y, *z;
429
430 n1 = npt1->len;
431 if(n1 == 0)
432 return(0);
433 n2 = npt2->len;
434 if(n2 == 0)
435 return(aERROR); /* Division durch 0 */
436
437 x = hilf + prec + 2;
438 y = x + (prec << 1);
439 z = y + prec;
440
441 diff = n2 - prec;
442 if(diff > 0) {
443 y = npt2->digits + diff;
444 npt2->expo += diff << 4;
445 n2 -= diff;
446 }
447 else
448 y = npt2->digits;
449
450 diff = n2 + prec - n1;
451 if(diff > 0) {
452 setarr(x,(int)diff,0);
453 cpyarr(npt1->digits,n1,x+diff);
454 }
455 else
456 cpyarr(npt1->digits-diff,(int)(n1+diff),x);
457 n1 += diff;
458 npt1->expo -= diff << 4;
459 n = divbig(x,n1,y,n2,z,&rlen,hilf);
460 diff = n - prec; /* diff is 0 or 1 */
461 n -= diff;
462 cpyarr(z+diff,n,npt1->digits);
463 exponent = npt1->expo - npt2->expo + (diff << 4);
464 if(exponent >= maxfltex)
465 return(aERROR);
466 else if(exponent <= -maxfltex) {
467 exponent = MOSTNEGEX;
468 n = 0;
469 }
470 npt1->expo = exponent;
471 npt1->len = n;
472 if(n == 0)
473 npt1->sign = 0;
474 else if(npt1->sign == npt2->sign || (npt1->sign && npt2->sign))
475 npt1->sign = 0;
476 else
477 npt1->sign = MINUSBYTE;
478 return(n);
479 }
480 /*-------------------------------------------------------------------*/
481 /*
482 ** die durch *nptr gegebene Zahl ist etwa 10**decdigs(nptr)
483 */
decdigs(nptr)484 PRIVATE long decdigs(nptr)
485 numdata *nptr;
486 {
487 static unsigned a[3] = {3,31,22088};
488 /* log(2)/log(10) = 1/3 - 1/31 - 1/22088 */
489 word4 u;
490 long b;
491 int n, sign;
492
493 b = n = nptr->len - 1;
494 if(n < 0)
495 return(MOSTNEGEX);
496 b <<= 4;
497 b += nptr->expo + bitlen(nptr->digits[n]);
498 sign = (b < 0);
499 u = (sign ? -b : b);
500 b = (u/a[0] - u/a[1] - u/a[2]);
501 b = (sign ? -b : b);
502 return(b);
503 }
504 /*-------------------------------------------------------------------*/
505 /*
506 ** berechnet base**a in *nptr mit prec 16-bit-Stellen
507 */
pwrtrunc(prec,base,a,nptr,hilf)508 PUBLIC int pwrtrunc(prec,base,a,nptr,hilf)
509 int prec;
510 unsigned base, a;
511 numdata *nptr;
512 word2 *hilf;
513 {
514 word2 *pow, *temp;
515 long offs, offs1;
516 unsigned bb;
517 int len;
518
519 temp = hilf + prec + 2;
520
521 nptr->sign = 0;
522 pow = nptr->digits;
523 pow[0] = (a ? base : 1);
524 len = 1;
525 offs = 0;
526 if(a <= 1) {
527 bb = 0;
528 }
529 else if(a <= 0xFFFF) {
530 bb = 1;
531 bb <<= (bitlen(a)-2);
532 }
533 else {
534 bb = 0x8000;
535 bb <<= (bitlen(a>>16)-1);
536 }
537 while(bb) {
538 len = multbig(pow,len,pow,len,temp,hilf);
539 offs += offs;
540 if((offs1 = len-prec) > 0) {
541 offs += offs1;
542 cpyarr(temp+offs1,prec,temp);
543 len = prec;
544 }
545 if(a & bb) {
546 len = multarr(temp,len,base,pow);
547 if(len > prec) { /* dann len = prec+1 */
548 offs++;
549 cpyarr(pow+1,prec,pow);
550 len = prec;
551 }
552 }
553 else
554 cpyarr(temp,len,pow);
555 bb >>= 1;
556 }
557 nptr->len = len;
558 nptr->expo = (offs << 4);
559 return(len);
560 }
561 /*-------------------------------------------------------------------*/
562 /*
563 ** Verwandelt float-Zahl (gegeben in *p) in bcd-Zahl mit places Stellen
564 ** Ist Rueckgabewert len != 0, so ist len = places und das
565 ** Ergebnis gleich
566 ** (arr,len) * 10**offs
567 ** wobei arr = nptr->digits, offs = nptr->expo.
568 ** Rueckgabewert 0 bedeutet die Zahl 0.0
569 ** arr muss mindestens die word2-Laenge (2 + places/4) haben
570 ** Platz hilf muss word2-Laenge >= places + places/4 + 8 haben
571 */
float2bcd(places,p,nptr,hilf)572 PUBLIC int float2bcd(places,p,nptr,hilf)
573 int places;
574 truc *p;
575 numdata *nptr;
576 word2 *hilf;
577 {
578 word2 *arr;
579 long d, d1;
580 int prec, prec1, n, sh;
581 int flg = *FLAGPTR(p);
582
583 prec = fltprec(flg);
584 prec1 = places/4;
585 if(prec > prec1)
586 prec = prec1;
587 prec +=1;
588
589 if(getnumtrunc(prec,p,nptr) == 0)
590 return(0);
591 d = decdigs(nptr);
592 d1 = places - d + 3; /* +3 wegen Rundungsfehlern */
593
594 n = malzehnhoch(prec,nptr,d1,hilf);
595 arr = nptr->digits;
596 cpyarr(arr,n,hilf);
597 sh = nptr->expo; /* sh ist int! */
598 n = shiftarr(hilf,n,sh);
599 nptr->len = big2bcd(hilf,n,arr);
600 nptr->expo = -d1;
601 return(roundbcd(places,nptr));
602 }
603 /*--------------------------------------------------------------*/
604 /*
605 ** Die in *nptr gegebene Float-Zahl wird auf prec signifikante
606 ** Stellen gerundet
607 */
roundbcd(prec,nptr)608 PUBLIC int roundbcd(prec,nptr)
609 int prec;
610 numdata *nptr;
611 {
612 word2 *arr;
613 int sh, len;
614 int carry, dig;
615
616 len = nptr->len;
617 sh = len - prec;
618 if(sh <= 0 || prec <= 0)
619 return(len);
620 arr = nptr->digits;
621 dig = nibdigit(arr,sh-1);
622 carry = (dig >= 5 ? 1 : 0);
623 len = shiftbcd(arr,len,-sh);
624 nptr->expo += sh;
625 if(carry)
626 len = incbcd(arr,len,carry);
627 if(len > prec) {
628 len = shiftbcd(arr,len,-1);
629 nptr->expo++;
630 }
631 return(nptr->len = len);
632 }
633 /*------------------------------------------------------------------*/
634 /*
635 ** Multipliziert die in *nptr gegebene Zahl mit 10**d
636 ** auf prec 16-Bit-Stellen genau
637 */
malzehnhoch(prec,nptr,d,hilf)638 PRIVATE int malzehnhoch(prec,nptr,d,hilf)
639 int prec;
640 long d;
641 numdata *nptr;
642 word2 *hilf;
643 {
644 numdata p10;
645 word2 *hilf1;
646 unsigned int dd;
647 int n;
648
649 p10.digits = hilf;
650 hilf1 = hilf + prec + 2;
651
652 dd = (d >= 0 ? d : -d);
653 pwrtrunc(prec,10,dd,&p10,hilf1);
654 if(d >= 0)
655 n = multtrunc(prec,nptr,&p10,hilf1);
656 else
657 n = divtrunc(prec,nptr,&p10,hilf1);
658 return(n);
659 }
660 /*------------------------------------------------------------------*/
661 /*
662 ** verwandelt die in *nptr gegebene Zahl, wobei sich nptr->expo
663 ** auf die Basis 10 bezieht, in Zahl bezueglich Basis 2
664 */
flodec2bin(prec,nptr,hilf)665 PUBLIC int flodec2bin(prec,nptr,hilf)
666 int prec;
667 numdata *nptr;
668 word2 *hilf;
669 {
670 long d;
671
672 d = nptr->expo;
673 nptr->expo = 0;
674 return(malzehnhoch(prec,nptr,d,hilf));
675 }
676 /*------------------------------------------------------------------*/
677 /*
678 ** verwandelt integer x (darf nur 16 bit enthalten) in numdata
679 */
int2numdat(x,nptr)680 PUBLIC void int2numdat(x,nptr)
681 int x;
682 numdata *nptr;
683 {
684 if(x < 0) {
685 x = -x;
686 nptr->sign = MINUSBYTE;
687 }
688 else
689 nptr->sign = 0;
690 *nptr->digits = x;
691 nptr->len = (x ? 1 : 0);
692 nptr->expo = (x ? 0 : MOSTNEGEX);
693 }
694 /*------------------------------------------------------------------*/
cpynumdat(npt1,npt2)695 PUBLIC void cpynumdat(npt1,npt2)
696 numdata *npt1, *npt2;
697 {
698 npt2->sign = npt1->sign;
699 npt2->len = npt1->len;
700 cpyarr(npt1->digits,npt1->len,npt2->digits);
701 npt2->expo = npt1->expo;
702 }
703 /*-------------------------------------------------------------------*/
704 /*
705 ** Ergibt +1, 0, -1, je nachdem Zahl in *ptr positiv, null oder
706 ** negativ ist.
707 */
numposneg(ptr)708 PUBLIC int numposneg(ptr)
709 truc *ptr;
710 {
711 int flg = *FLAGPTR(ptr);
712
713 if(((flg == fFIXNUM) && *SIGNPTR(ptr)) ||
714 ((flg == fBIGNUM) && *SIGNUMPTR(ptr)) ||
715 ((flg >= fFLTOBJ) && (*SIGNUMPTR(ptr) & FSIGNBIT)))
716 return(-1);
717 else if((*ptr == zero) || ((flg >= fFLTOBJ) && (flg & FLTZEROBIT)))
718 return(0);
719 else
720 return(1);
721 }
722 /*------------------------------------------------------------------*/
723 /*
724 ** Loescht destruktiv das Vorzeichen der in *ptr gegebenen Zahl.
725 */
wipesign(ptr)726 PUBLIC truc wipesign(ptr)
727 truc *ptr;
728 {
729 int flg = *FLAGPTR(ptr);
730
731 if(flg == fFIXNUM)
732 *SIGNPTR(ptr) = 0;
733 else if(flg >= fFLTOBJ) {
734 if((flg & FLTZEROBIT) == 0)
735 *SIGNUMPTR(ptr) &= ~FSIGNBIT;
736 }
737 else /* bignums */
738 *SIGNUMPTR(ptr) = 0;
739 return(*ptr);
740 }
741 /*------------------------------------------------------------------*/
742 /*
743 ** Aendert destruktiv das Vorzeichen der in *ptr gegebenen Zahl.
744 */
changesign(ptr)745 PUBLIC truc changesign(ptr)
746 truc *ptr;
747 {
748 int sign, flg;
749
750 sign = numposneg(ptr);
751 if(sign == 0)
752 return *ptr;
753 else if(sign < 0)
754 return wipesign(ptr);
755
756 /* now *ptr is positive */
757 flg = *FLAGPTR(ptr);
758 if(flg == fFIXNUM)
759 *SIGNPTR(ptr) = MINUSBYTE;
760 else if(flg >= fFLTOBJ) {
761 *SIGNUMPTR(ptr) |= FSIGNBIT;
762 }
763 else /* bignums */
764 *SIGNUMPTR(ptr) = MINUSBYTE;
765 return(*ptr);
766 }
767 /*------------------------------------------------------------------*/
768 /*
769 ** holt long aus *ptr;
770 ** falls abs(Zahl) >= 2 ** 31, Returnwert: LONGERROR
771 */
intretr(ptr)772 PUBLIC long intretr(ptr)
773 truc *ptr;
774 {
775 word2 *x;
776 word4 u;
777 long res;
778 int n, sign;
779
780 n = bigref(ptr,&x,&sign);
781 if(n < 0 || n > 2 || (u = big2long(x,n)) >= 0x80000000)
782 return(LONGERROR);
783 res = u;
784 return(sign ? -res : res);
785 }
786 /*------------------------------------------------------------------*/
787 /*
788 ** Erzeugt Referenz auf Integer (bignum, fixnum oder gf2n_int) ohne Kopie
789 ** Rueckgabe aERROR, falls *ptr kein fixnum, bignum oder gf2nint
790 */
bigref(ptr,xp,sp)791 PUBLIC int bigref(ptr,xp,sp)
792 truc *ptr;
793 word2 **xp;
794 int *sp;
795 {
796 struct bigcell *big;
797 word2 *x;
798 int flg = *FLAGPTR(ptr);
799
800 if(flg == fFIXNUM) {
801 *xp = x = WORD2PTR(ptr);
802 *sp = *SIGNPTR(ptr);
803 return(*x ? 1 : 0);
804 }
805 else if(flg == fBIGNUM || flg == fGF2NINT) {
806 big = (struct bigcell *)TAddress(ptr);
807 *sp = big->signum;
808 *xp = &(big->digi0);
809 return(big->len);
810 }
811 else
812 return(aERROR);
813 }
814 /*------------------------------------------------------------------*/
815 /*
816 ** holt Integer (fixnum oder bignum) aus *ptr nach x
817 ** Vorzeichen in *sp
818 ** Rueckgabe aERROR, falls *ptr kein fixnum oder bignum
819 */
bigretr(ptr,x,sp)820 PUBLIC int bigretr(ptr,x,sp)
821 truc *ptr;
822 word2 *x;
823 int *sp;
824 {
825 word2 *z;
826 int n;
827
828 n = bigref(ptr,&z,sp);
829 if(n != aERROR)
830 cpyarr(z,n,x);
831 return(n);
832 }
833 /*-------------------------------------------------------------*/
834 /*
835 ** holt Integer (fixnum oder bignum) aus *ptr nach x
836 ** in Zweier-Komplement-Darstellung
837 ** x[n] = 0, falls positiv;
838 ** x[n] = 0xFFFF, falls negativ.
839 */
twocretr(ptr,x)840 PUBLIC int twocretr(ptr,x)
841 truc *ptr;
842 word2 *x;
843 {
844 int sign, n;
845
846 n = bigretr(ptr,x,&sign);
847 if(n < 0)
848 return(n);
849 if(sign) {
850 n = decarr(x,n,1);
851 notarr(x,n);
852 x[n] = 0xFFFF;
853 }
854 else
855 x[n] = 0;
856 return(n);
857 }
858 /*-------------------------------------------------------------*/
twocadjust(x,n)859 PRIVATE int twocadjust(x,n)
860 word2 *x;
861 int n;
862 {
863 word2 u = x[n];
864
865 while((n > 0) && (x[n-1] == u))
866 n--;
867 return(n);
868 }
869 /*-------------------------------------------------------------*/
870 /*
871 ** Bitwise and of two bignums (x,n), (y,m)
872 ** in two's complement representation
873 ** Destructively replaces (x,n) by result
874 */
and2arr(x,n,y,m)875 PUBLIC int and2arr(x,n,y,m)
876 word2 *x,*y;
877 int n,m;
878 {
879 if(n < m && x[n]) {
880 setarr(x+n+1,m-n,0xFFFF);
881 n = m;
882 }
883 else if(m < n && !y[m]) {
884 n = m;
885 }
886 andarr(x,m+1,y);
887 n = twocadjust(x,n);
888 return(n);
889 }
890 /*-------------------------------------------------------------*/
891 /*
892 ** Bitwise or of two bignums (x,n), (y,m)
893 ** in two's complement representation
894 ** Destructively replaces (x,n) by result
895 */
or2arr(x,n,y,m)896 PUBLIC int or2arr(x,n,y,m)
897 word2 *x,*y;
898 int n,m;
899 {
900 if(n < m && !x[n]) {
901 setarr(x+n+1,m-n,0);
902 n = m;
903 }
904 else if(m < n && y[m]) {
905 n = m;
906 }
907 orarr(x,m+1,y);
908 n = twocadjust(x,n);
909 return(n);
910 }
911 /*-------------------------------------------------------------*/
912 /*
913 ** Bitwise xor of two bignums (x,n), (y,m)
914 ** in two's complement representation
915 ** Destructively replaces (x,n) by result
916 */
xor2arr(x,n,y,m)917 PUBLIC int xor2arr(x,n,y,m)
918 word2 *x,*y;
919 int n,m;
920 {
921 if(n < m) {
922 setarr(x+n+1,m-n,x[n]);
923 n = m;
924 }
925 else if(m < n && y[m]) {
926 setarr(y+m+1,n-m,0xFFFF);
927 m = n;
928 }
929 xorarr(x,m+1,y);
930 n = twocadjust(x,n);
931 return(n);
932 }
933 /*-------------------------------------------------------------*/
934 /*
935 ** Bitwise xor of two bitvectors (x,n), (y,m)
936 ** Destructively replaces (x,n) by result
937 */
xorbitvec(x,n,y,m)938 PUBLIC int xorbitvec(x,n,y,m)
939 word2 *x,*y;
940 int n,m;
941 {
942 if(n < m) {
943 setarr(x+n,m-n,0);
944 n = m;
945 }
946 xorarr(x,m,y);
947 while(n > 0 && x[n-1] == 0)
948 n--;
949 return(n);
950 }
951 /*-------------------------------------------------------------*/
bit_length(x,n)952 PUBLIC long bit_length(x,n)
953 word2 *x;
954 int n;
955 {
956 long b;
957
958 if(n > 0) {
959 b = n - 1;
960 b <<= 4;
961 b += bitlen(x[n-1]);
962 }
963 else
964 b = 0;
965 return(b);
966 }
967 /*-------------------------------------------------------------*/
968 /*
969 ** Prueft, ob die Elemente argptr[i], 0 <= i < n, Integers sind
970 ** Rueckgabewert: fFIXNUM, fBIGNUM oder aERROR
971 */
chkints(sym,argptr,n)972 PUBLIC int chkints(sym,argptr,n)
973 truc sym;
974 truc *argptr;
975 int n;
976 {
977 int flg, flg0 = fFIXNUM;
978
979 while(--n >= 0) {
980 flg = *FLAGPTR(argptr);
981 if(flg < fFIXNUM || flg > fBIGNUM)
982 return(error(sym,err_int,*argptr));
983 else if(flg > flg0)
984 flg0 = flg;
985 argptr++;
986 }
987 return(flg0);
988 }
989 /*------------------------------------------------------------------*/
990 /*
991 ** checks whether element *ptr is an integer
992 ** Return value: fFIXNUM, fBIGNUM or aERROR
993 */
chkint(sym,ptr)994 PUBLIC int chkint(sym,ptr)
995 truc sym;
996 truc *ptr;
997 {
998 int flg;
999
1000 flg = *FLAGPTR(ptr);
1001 if(flg < fFIXNUM || flg > fBIGNUM)
1002 return(error(sym,err_int,*ptr));
1003 else
1004 return flg;
1005 }
1006 /*------------------------------------------------------------------*/
1007 /*
1008 ** checks whether element *ptr is an integer or gf2nint
1009 ** Return value: fFIXNUM, fBIGNUM, fGF2NINT or aERROR
1010 */
chkintt(sym,ptr)1011 PUBLIC int chkintt(sym,ptr)
1012 truc sym;
1013 truc *ptr;
1014 {
1015 int flg;
1016
1017 flg = *FLAGPTR(ptr);
1018 if(flg < fINTTYPE0 || flg > fINTTYPE1)
1019 return(error(sym,err_intt,*ptr));
1020 else
1021 return flg;
1022 }
1023 /*------------------------------------------------------------------*/
1024 /*
1025 ** checks whether element *ptr is a nonzero integer
1026 ** Return value: fFIXNUM, fBIGNUM or aERROR
1027 */
chkintnz(sym,ptr)1028 PUBLIC int chkintnz(sym,ptr)
1029 truc sym;
1030 truc *ptr;
1031 {
1032 int flg;
1033
1034 flg = *FLAGPTR(ptr);
1035 if(flg < fFIXNUM || flg > fBIGNUM)
1036 return(error(sym,err_int,*ptr));
1037 else if(*ptr == zero)
1038 return(error(sym,err_div,voidsym));
1039 else
1040 return flg;
1041 }
1042 /*-------------------------------------------------------------*/
1043 /*
1044 ** Argument *vptr must be a vector.
1045 ** The function checks whether all components of this
1046 ** vector are integers.
1047 ** Return value: fFIXNUM, fBIGNUM or aERROR
1048 */
chkintvec(sym,vptr)1049 PUBLIC int chkintvec(sym,vptr)
1050 truc sym;
1051 truc *vptr;
1052 {
1053 struct vector *vec;
1054 truc *ptr;
1055 int len;
1056 int flg, flg0 = fFIXNUM;
1057
1058 vec = (struct vector *)TAddress(vptr);
1059 len = vec->len;
1060 ptr = &(vec->ele0);
1061
1062 while(--len >= 0) {
1063 flg = *FLAGPTR(ptr);
1064 if(flg < fFIXNUM || flg > fBIGNUM)
1065 return(error(sym,err_int,*ptr));
1066 else if(flg > flg0)
1067 flg0 = flg;
1068 ptr++;
1069 }
1070 return(flg0);
1071 }
1072 /*-------------------------------------------------------------*/
1073 /*
1074 ** Prueft, ob die Elemente argptr[i], 0 <= i < n, Integers oder
1075 ** Floats sind
1076 ** Rueckgabewert: fFIXNUM, fBIGNUM, Float-Flag oder aERROR
1077 */
chknums(sym,argptr,n)1078 PUBLIC int chknums(sym,argptr,n)
1079 truc sym;
1080 truc *argptr;
1081 int n;
1082 {
1083 int flg, flg0 = fFIXNUM;
1084
1085 while(--n >= 0) {
1086 flg = *FLAGPTR(argptr);
1087 if(flg < fFIXNUM)
1088 return(error(sym,err_num,*argptr));
1089 else if(flg > flg0)
1090 flg0 = flg;
1091 argptr++;
1092 }
1093 return(flg0);
1094 }
1095 /*-------------------------------------------------------------*/
chknum(sym,ptr)1096 PUBLIC int chknum(sym,ptr)
1097 truc sym;
1098 truc *ptr;
1099 {
1100 int flg = *FLAGPTR(ptr);
1101
1102 if(flg < fFIXNUM)
1103 return(error(sym,err_num,*ptr));
1104 return(flg);
1105 }
1106 /*-------------------------------------------------------------*/
1107 /*
1108 ** Argument *vptr must be a vector.
1109 ** The function checks whether all components of this
1110 ** vector are integers.
1111 ** Return value: fFIXNUM, fBIGNUM, Float-flag or aERROR
1112 */
chknumvec(sym,vptr)1113 PUBLIC int chknumvec(sym,vptr)
1114 truc sym;
1115 truc *vptr;
1116 {
1117 struct vector *vec;
1118 truc *ptr;
1119 int len;
1120 int flg, flg0 = fFIXNUM;
1121
1122 vec = (struct vector *)TAddress(vptr);
1123 len = vec->len;
1124 ptr = &(vec->ele0);
1125
1126 while(--len >= 0) {
1127 flg = *FLAGPTR(ptr);
1128 if(flg < fFIXNUM)
1129 return(error(sym,err_num,*ptr));
1130 else if(flg > flg0)
1131 flg0 = flg;
1132 ptr++;
1133 }
1134 return(flg0);
1135 }
1136 /***************************************************************/
1137