]>
Commit | Line | Data |
---|---|---|
0095f096 | 1 | /* Copyright Digital Equipment Corporation & INRIA 1988, 1989 */ |
2 | /* Last modified_on Fri Mar 2 16:49:23 GMT+1:00 1990 by herve */ | |
3 | /* modified_on Mon Feb 19 20:22:20 GMT+1:00 1990 by shand */ | |
4 | ||
5 | ||
6 | /* KerN.c: the kernel written in C */ | |
7 | ||
8 | /* | |
9 | * Description of types and constants. | |
10 | * | |
11 | * Several conventions are used in the commentary: | |
12 | * A "BigNum" is the name for an infinite-precision number. | |
13 | * Capital letters (e.g., "N") are used to refer to the value of BigNums. | |
14 | * The word "digit" refers to a single BigNum digit. | |
15 | * The notation "Size(N)" refers to the number of digits in N, | |
16 | * which is typically passed to the subroutine as "nl". | |
17 | * The notation "Length(N)" refers to the number of digits in N, | |
18 | * not including any leading zeros. | |
19 | * The word "Base" is used for the number 2 ** BN_DIGIT_SIZE, where | |
20 | * BN_DIGIT_SIZE is the number of bits in a single BigNum digit. | |
21 | * The expression "BBase(N)" is used for Base ** NumDigits(N). | |
22 | * The term "leading zeros" refers to any zeros before the most | |
23 | * significant digit of a number. | |
24 | * | |
25 | * | |
26 | * In the code, we have: | |
27 | * | |
28 | * "nn" is a pointer to a big number, | |
29 | * "nl" is the number of digits from nn, | |
30 | * "d" is a digit. | |
31 | * | |
32 | */ | |
33 | ||
34 | ||
35 | /*\f*/ | |
36 | ||
37 | ||
38 | #define BNNMACROS_OFF | |
39 | #include "BigNum.h" | |
40 | ||
41 | ||
42 | /*** copyright ***/ | |
43 | ||
44 | static char copyright[]="@(#)KerN.c: copyright Digital Equipment Corporation & INRIA 1988, 1989\n"; | |
45 | ||
46 | ||
47 | /******* non arithmetic access to digits ********/ | |
48 | ||
49 | ||
50 | void BnnSetToZero (nn, nl) | |
51 | ||
52 | register BigNum nn; | |
53 | register int nl; | |
54 | ||
55 | /* | |
56 | * Sets all the specified digits of the BigNum to 0 | |
57 | */ | |
58 | ||
59 | { | |
60 | #ifdef NOMEM | |
61 | while (--nl >= 0) | |
62 | *(nn++) = 0; | |
63 | #else | |
64 | memset (nn, 0, nl*BN_DIGIT_SIZE/BN_BYTE_SIZE); | |
65 | #endif | |
66 | } | |
67 | ||
68 | /***************************************/ | |
69 | ||
70 | ||
71 | void BnnAssign (mm, nn, nl) | |
72 | ||
73 | register BigNum mm, nn; | |
74 | register int nl; | |
75 | ||
76 | /* | |
77 | * Copies N => M | |
78 | */ | |
79 | ||
80 | { | |
81 | if (mm < nn || mm > nn+nl) | |
82 | #ifdef NOMEM | |
83 | while (--nl >= 0) | |
84 | *mm++ = *nn++; | |
85 | #else | |
86 | /* be care: bcopy (SRC, DEST, L): SRC-->DEST !!! */ | |
0095f096 | 87 | memmove (mm, nn, nl*BN_DIGIT_SIZE/BN_BYTE_SIZE); |
0095f096 | 88 | #endif |
89 | else | |
90 | if (mm > nn) | |
91 | { | |
92 | nn += nl; | |
93 | mm += nl; | |
94 | while (--nl >= 0) | |
95 | *--mm = *--nn; | |
96 | } | |
97 | } | |
98 | ||
99 | /***************************************/ | |
100 | /*\f*/ | |
101 | ||
102 | ||
103 | void BnnSetDigit (nn, d) | |
104 | ||
105 | BigNum nn; | |
106 | int d; | |
107 | ||
108 | /* | |
109 | * Sets a single digit of N to the passed value | |
110 | */ | |
111 | ||
112 | { | |
113 | *nn = d; | |
114 | } | |
115 | ||
116 | /***************************************/ | |
117 | ||
118 | ||
119 | BigNumDigit BnnGetDigit (nn) | |
120 | ||
121 | BigNum nn; | |
122 | ||
123 | /* | |
124 | * Returns the single digit pointed by N | |
125 | */ | |
126 | ||
127 | { | |
128 | return (*nn); | |
129 | } | |
130 | ||
131 | /***************************************/ | |
132 | /*\f*/ | |
133 | ||
134 | ||
135 | BigNumLength BnnNumDigits (nn, nl) | |
136 | ||
137 | register BigNum nn; | |
138 | register int nl; | |
139 | ||
140 | /* | |
141 | * Returns the number of digits of N, not counting leading zeros | |
142 | */ | |
143 | ||
144 | { | |
145 | nn += nl; | |
146 | ||
147 | while (nl != 0 && *--nn == 0) | |
148 | nl--; | |
149 | ||
150 | return (nl == 0 ? 1 : nl); | |
151 | } | |
152 | ||
153 | /***************************************/ | |
154 | ||
155 | ||
156 | BigNumDigit BnnNumLeadingZeroBitsInDigit (d) | |
157 | ||
158 | BigNumDigit d; | |
159 | ||
160 | /* | |
161 | * Returns the number of leading zero bits in a digit | |
162 | */ | |
163 | ||
164 | { | |
165 | register BigNumDigit mask = 1 << (BN_DIGIT_SIZE - 1); | |
166 | register int p = 0; | |
167 | ||
168 | ||
169 | if (d == 0) | |
170 | return (BN_DIGIT_SIZE); | |
171 | ||
172 | while ((d & mask) == 0) | |
173 | { | |
174 | p++; | |
175 | mask >>= 1; | |
176 | } | |
177 | ||
178 | return (p); | |
179 | } | |
180 | ||
181 | /***************************************/ | |
182 | /*\f*/ | |
183 | ||
184 | /************** Predicates on one digit ***************/ | |
185 | ||
186 | ||
187 | Boolean BnnDoesDigitFitInWord (d) | |
188 | ||
189 | BigNumDigit d; | |
190 | ||
191 | /* | |
192 | * Returns TRUE iff the digit can be represented in just BN_WORD_SIZE bits | |
193 | */ | |
194 | { | |
195 | /* The C compiler must evaluate the predicate at compile time */ | |
196 | if (BN_DIGIT_SIZE > BN_WORD_SIZE) | |
197 | return (d >= 1 << BN_WORD_SIZE ? FALSE : TRUE); | |
198 | else | |
199 | return (TRUE); | |
200 | } | |
201 | ||
202 | /***************************************/ | |
203 | ||
204 | ||
205 | Boolean BnnIsDigitZero (d) | |
206 | ||
207 | BigNumDigit d; | |
208 | ||
209 | /* Returns TRUE iff digit = 0 */ | |
210 | ||
211 | { | |
212 | return (d == 0); | |
213 | } | |
214 | ||
215 | /***************************************/ | |
216 | ||
217 | ||
218 | Boolean BnnIsDigitNormalized (d) | |
219 | ||
220 | BigNumDigit d; | |
221 | ||
222 | /* | |
223 | * Returns TRUE iff Base/2 <= digit < Base | |
224 | * i.e., if digit's leading bit is 1 | |
225 | */ | |
226 | ||
227 | { | |
228 | return (d & (1 << (BN_DIGIT_SIZE - 1)) ? TRUE : FALSE); | |
229 | } | |
230 | ||
231 | /***************************************/ | |
232 | ||
233 | ||
234 | Boolean BnnIsDigitOdd (d) | |
235 | ||
236 | BigNumDigit d; | |
237 | ||
238 | /* | |
239 | * Returns TRUE iff digit is odd | |
240 | */ | |
241 | ||
242 | { | |
243 | return (d & 1 ? TRUE : FALSE); | |
244 | } | |
245 | ||
246 | /***************************************/ | |
247 | ||
248 | ||
249 | BigNumCmp BnnCompareDigits (d1, d2) | |
250 | ||
251 | BigNumDigit d1, d2; | |
252 | ||
253 | /* | |
254 | * Returns BN_GREATER if digit1 > digit2 | |
255 | * BN_EQUAL if digit1 = digit2 | |
256 | * BN_LESS if digit1 < digit2 | |
257 | */ | |
258 | ||
259 | { | |
260 | return (d1 > d2 ? BN_GT : (d1 == d2 ? BN_EQ : BN_LT)); | |
261 | } | |
262 | ||
263 | /***************** Logical operations ********************/ | |
264 | ||
265 | ||
266 | void BnnComplement (nn, nl) | |
267 | ||
268 | register BigNum nn; | |
269 | register int nl; | |
270 | ||
271 | /* | |
272 | * Performs the computation BBase(N) - N - 1 => N | |
273 | */ | |
274 | ||
275 | { | |
276 | while (--nl >= 0) | |
277 | *(nn++) ^= -1; | |
278 | } | |
279 | ||
280 | /***************************************/ | |
281 | /*\f*/ | |
282 | ||
283 | ||
284 | void BnnAndDigits (n, d) | |
285 | ||
286 | BigNum n; | |
287 | BigNumDigit d; | |
288 | ||
289 | /* | |
290 | * Returns the logical computation n[0] AND d in n[0] | |
291 | */ | |
292 | ||
293 | { | |
294 | *n &= d; | |
295 | } | |
296 | ||
297 | /***************************************/ | |
298 | ||
299 | ||
300 | void BnnOrDigits (n, d) | |
301 | ||
302 | BigNum n; | |
303 | BigNumDigit d; | |
304 | ||
305 | /* | |
306 | * Returns the logical computation n[0] OR d2 in n[0]. | |
307 | */ | |
308 | ||
309 | { | |
310 | *n |= d; | |
311 | } | |
312 | ||
313 | /***************************************/ | |
314 | ||
315 | ||
316 | void BnnXorDigits (n, d) | |
317 | ||
318 | BigNum n; | |
319 | BigNumDigit d; | |
320 | ||
321 | /* | |
322 | * Returns the logical computation n[0] XOR d in n[0]. | |
323 | */ | |
324 | ||
325 | { | |
326 | *n ^= d; | |
327 | } | |
328 | ||
329 | /***************************************/ | |
330 | /*\f*/ | |
331 | ||
332 | /****************** Shift operations *******************/ | |
333 | ||
334 | ||
335 | BigNumDigit BnnShiftLeft (mm, ml, nbits) | |
336 | ||
337 | register BigNum mm; | |
338 | register int ml; | |
339 | int nbits; | |
340 | ||
341 | /* | |
342 | * Shifts M left by "nbits", filling with 0s. | |
343 | * Returns the leftmost "nbits" of M in a digit. | |
344 | * Assumes 0 <= nbits < BN_DIGIT_SIZE. | |
345 | */ | |
346 | ||
347 | { | |
348 | register BigNumDigit res = 0, save; | |
349 | int rnbits; | |
350 | ||
351 | ||
352 | if (nbits != 0) | |
353 | { | |
354 | rnbits = BN_DIGIT_SIZE - nbits; | |
355 | ||
356 | while (--ml >= 0) | |
357 | { | |
358 | save = *mm; | |
359 | *mm++ = (save << nbits) | res; | |
360 | res = save >> rnbits; | |
361 | } | |
362 | } | |
363 | ||
364 | return (res); | |
365 | } | |
366 | ||
367 | /***************************************/ | |
368 | /*\f*/ | |
369 | ||
370 | ||
371 | BigNumDigit BnnShiftRight (mm, ml, nbits) | |
372 | ||
373 | register BigNum mm; | |
374 | register int ml; | |
375 | int nbits; | |
376 | ||
377 | /* | |
378 | * Shifts M right by "nbits", filling with 0s. | |
379 | * Returns the rightmost "nbits" of M in a digit. | |
380 | * Assumes 0 <= nbits < BN_DIGIT_SIZE. | |
381 | */ | |
382 | ||
383 | { | |
384 | register BigNumDigit res = 0, save; | |
385 | int lnbits; | |
386 | ||
387 | ||
388 | if (nbits != 0) | |
389 | { | |
390 | mm += ml; | |
391 | lnbits = BN_DIGIT_SIZE - nbits; | |
392 | ||
393 | while (--ml >= 0) | |
394 | { | |
395 | save = *(--mm); | |
396 | *mm = (save >> nbits) | res; | |
397 | res = save << lnbits; | |
398 | } | |
399 | } | |
400 | ||
401 | return (res); | |
402 | } | |
403 | ||
404 | /***************************************/ | |
405 | /*\f*/ | |
406 | ||
407 | ||
408 | /******************* Additions **************************/ | |
409 | ||
410 | ||
411 | BigNumCarry BnnAddCarry (nn, nl, carryin) | |
412 | ||
413 | register BigNum nn; | |
414 | register int nl; | |
415 | BigNumCarry carryin; | |
416 | ||
417 | /* | |
418 | * Performs the sum N + CarryIn => N. | |
419 | * Returns the CarryOut. | |
420 | */ | |
421 | ||
422 | { | |
423 | if (carryin == 0) | |
424 | return (0); | |
425 | ||
426 | if (nl == 0) | |
427 | return (1); | |
428 | ||
429 | while (--nl >= 0 && !(++(*nn++))) | |
430 | ; | |
431 | ||
432 | return (nl >= 0 ? 0 : 1); | |
433 | } | |
434 | ||
435 | /***************************************/ | |
436 | /*\f*/ | |
437 | ||
438 | ||
439 | BigNumCarry BnnAdd (mm, ml, nn, nl, carryin) | |
440 | ||
441 | register BigNum mm, nn; | |
442 | int ml; | |
443 | register int nl; | |
444 | BigNumCarry carryin; | |
445 | ||
446 | /* | |
447 | * Performs the sum M + N + CarryIn => M. | |
448 | * Returns the CarryOut. | |
449 | * Assumes Size(M) >= Size(N). | |
450 | */ | |
451 | ||
452 | { | |
453 | register BigNumProduct c = carryin; | |
454 | ||
455 | ||
456 | ml -= nl; | |
457 | /* test computed at compile time */ | |
458 | if (sizeof (BigNumProduct) > sizeof (BigNumDigit)) | |
459 | { | |
460 | while (--nl >= 0) | |
461 | { | |
462 | c += *mm + *(nn++); | |
463 | *(mm++) = c; | |
464 | c >>= BN_DIGIT_SIZE; | |
465 | } | |
466 | } | |
467 | else | |
468 | { | |
469 | register BigNumProduct save; | |
470 | ||
471 | while (--nl >= 0) | |
472 | { | |
473 | save = *mm; | |
474 | c += save; | |
475 | if (c < save) | |
476 | { | |
477 | *(mm++) = *(nn++); | |
478 | c = 1; | |
479 | } | |
480 | else | |
481 | { | |
482 | save = *(nn++); | |
483 | c += save; | |
484 | *(mm++) = c; | |
485 | c = (c < save) ? 1 : 0; | |
486 | } | |
487 | } | |
488 | } | |
489 | ||
490 | return (BnnAddCarry (mm, ml, (BigNumCarry) c)); | |
491 | } | |
492 | ||
493 | /***************************************/ | |
494 | /*\f*/ | |
495 | ||
496 | /****************** Subtraction *************************/ | |
497 | ||
498 | ||
499 | ||
500 | BigNumCarry BnnSubtractBorrow (nn, nl, carryin) | |
501 | ||
502 | register BigNum nn; | |
503 | register int nl; | |
504 | BigNumCarry carryin; | |
505 | ||
506 | /* | |
507 | * Performs the difference N + CarryIn - 1 => N. | |
508 | * Returns the CarryOut. | |
509 | */ | |
510 | ||
511 | { | |
512 | if (carryin == 1) | |
513 | return (1); | |
514 | if (nl == 0) | |
515 | return (0); | |
516 | ||
517 | #ifdef ibm032 | |
518 | while (--nl >= 0 && !(*nn)--) nn++; | |
519 | #else | |
520 | while (--nl >= 0 && !((*nn++)--)) | |
521 | #endif | |
522 | ; | |
523 | ||
524 | return (nl >= 0 ? 1 : 0); | |
525 | } | |
526 | ||
527 | /***************************************/ | |
528 | /*\f*/ | |
529 | ||
530 | ||
531 | BigNumCarry BnnSubtract (mm, ml, nn, nl, carryin) | |
532 | ||
533 | register BigNum mm, nn; | |
534 | int ml; | |
535 | register int nl; | |
536 | BigNumCarry carryin; | |
537 | ||
538 | /* | |
539 | * Performs the difference M - N + CarryIn - 1 => M. | |
540 | * Returns the CarryOut. | |
541 | * Assumes Size(M) >= Size(N). | |
542 | */ | |
543 | ||
544 | { | |
545 | register BigNumProduct c = carryin; | |
546 | register BigNumDigit invn; | |
547 | ||
548 | ||
549 | ml -= nl; | |
550 | /* test computed at compile time */ | |
551 | if (sizeof (BigNumProduct) > sizeof (BigNumDigit)) | |
552 | { | |
553 | while (--nl >= 0) | |
554 | { | |
555 | invn = *(nn++) ^ -1; | |
556 | c += *mm + invn; | |
557 | *(mm++) = c; | |
558 | c >>= BN_DIGIT_SIZE; | |
559 | } | |
560 | } | |
561 | else | |
562 | { | |
563 | register BigNumProduct save; | |
564 | ||
565 | while (--nl >= 0) | |
566 | { | |
567 | save = *mm; | |
568 | invn = *(nn++) ^ -1; | |
569 | c += save; | |
570 | ||
571 | if (c < save) | |
572 | { | |
573 | *(mm++) = invn; | |
574 | c = 1; | |
575 | } | |
576 | else | |
577 | { | |
578 | c += invn; | |
579 | *(mm++) = c; | |
580 | c = (c < invn) ? 1 : 0; | |
581 | } | |
582 | } | |
583 | } | |
584 | ||
585 | return (BnnSubtractBorrow (mm, ml, (BigNumCarry) c)); } | |
586 | ||
587 | ||
588 | /***************************************/ | |
589 | /*\f */ | |
590 | ||
591 | /***************** Multiplication ************************/ | |
592 | ||
593 | ||
594 | BigNumCarry BnnMultiplyDigit (pp, pl, mm, ml, d) | |
595 | ||
596 | register BigNum pp, mm; | |
597 | int pl, ml; | |
598 | BigNumDigit d; | |
599 | ||
600 | /* | |
601 | * Performs the product: | |
602 | * Q = P + M * d | |
603 | * BB = BBase(P) | |
604 | * Q mod BB => P | |
605 | * Q div BB => CarryOut | |
606 | * Returns the CarryOut. | |
607 | * Assumes Size(P) >= Size(M) + 1. | |
608 | */ | |
609 | ||
610 | { | |
611 | register BigNumProduct c = 0; | |
612 | ||
613 | ||
614 | if (d == 0) | |
615 | return (0); | |
616 | ||
617 | if (d == 1) | |
618 | return (BnnAdd (pp, pl, mm, ml, (BigNumCarry) 0)); | |
619 | ||
620 | pl -= ml; | |
621 | /* test computed at compile time */ | |
622 | if (sizeof (BigNumProduct) > sizeof (BigNumDigit)) | |
623 | { | |
624 | while (ml != 0) | |
625 | { | |
626 | ml--; | |
627 | c += *pp + (d * (*(mm++))); | |
628 | *(pp++) = c; | |
629 | c >>= BN_DIGIT_SIZE; | |
630 | } | |
631 | ||
632 | while (pl != 0) { | |
633 | pl--; | |
634 | c += *pp; | |
635 | *(pp++) = c; | |
636 | c >>= BN_DIGIT_SIZE; | |
637 | } | |
638 | ||
639 | return (c); | |
640 | } | |
641 | else | |
642 | { | |
643 | /* help for stupid compilers--may actually be counter | |
644 | productive on pipelined machines with decent register allocation!! */ | |
645 | #define m_digit X0 | |
646 | #define X3 Lm | |
647 | #define X1 Hm | |
648 | register BigNumDigit Lm, Hm, Ld, Hd, X0, X2 /*, X1, X3 */; | |
649 | ||
650 | Ld = d & ((1 << (BN_DIGIT_SIZE / 2)) -1); | |
651 | Hd = d >> (BN_DIGIT_SIZE / 2); | |
652 | while (ml != 0) | |
653 | { | |
654 | ml--; | |
655 | m_digit = *mm++; | |
656 | Lm = m_digit & ((1 << (BN_DIGIT_SIZE / 2)) -1); | |
657 | Hm = m_digit >> (BN_DIGIT_SIZE / 2); | |
658 | X0 = Ld * Lm; | |
659 | X2 = Hd * Lm; | |
660 | X3 = Hd * Hm; | |
661 | X1 = Ld * Hm; | |
662 | ||
663 | if ((c += X0) < X0) X3++; | |
664 | if ((X1 += X2) < X2) X3 += (1<<(BN_DIGIT_SIZE / 2)); | |
665 | X3 += (X1 >> (BN_DIGIT_SIZE / 2)); | |
666 | X1 <<= (BN_DIGIT_SIZE / 2); | |
667 | if ((c += X1) < X1) X3++; | |
668 | if ((*pp += c) < c) X3++; | |
669 | pp++; | |
670 | ||
671 | c = X3; | |
672 | #undef m_digit | |
673 | #undef X1 | |
674 | #undef X3 | |
675 | } | |
676 | ||
677 | X0 = *pp; | |
678 | c += X0; | |
679 | *(pp++) = c; | |
680 | ||
681 | if (c >= X0) | |
682 | return (0); | |
683 | ||
684 | pl--; | |
685 | while (pl != 0 && !(++(*pp++))) | |
686 | pl--; | |
687 | ||
688 | return (pl != 0 ? 0 : 1); | |
689 | } | |
690 | } | |
691 | ||
692 | #ifdef mips | |
693 | BigNumCarry BnnMultiply2Digit (pp, pl, mm, ml, d0, d1) | |
694 | ||
695 | register BigNum pp, mm; | |
696 | register int pl, ml; | |
697 | BigNumDigit d0, d1; | |
698 | ||
699 | /* | |
700 | * Provided for compatibility with mips assembler implementation. | |
701 | * Performs the product: | |
702 | * Q = P + M * d0_d1 | |
703 | * BB = BBase(P) | |
704 | * Q mod BB => P | |
705 | * Q div BB => CarryOut | |
706 | * Returns the CarryOut. | |
707 | * Assumes Size(P) >= Size(M) + 1. | |
708 | */ | |
709 | ||
710 | { | |
711 | return | |
712 | BnnMultiplyDigit (pp, pl, mm, ml, d0) | |
713 | + BnnMultiplyDigit (pp+1, pl-1, mm, ml, d1); | |
714 | } | |
715 | #endif /* mips */ | |
716 | ||
717 | /***************************************/ | |
718 | /*\f*/ | |
719 | ||
720 | /********************** Division *************************/ | |
721 | ||
722 | ||
723 | /* xh:xl -= yh:yl */ | |
724 | #define SUB(xh,xl,yh,yl) if (yl > xl) {xl -= yl; xh -= yh + 1;}\ | |
725 | else {xl -= yl; xh -= yh;} | |
726 | ||
727 | #define LOW(x) (x & ((1 << (BN_DIGIT_SIZE / 2)) -1)) | |
728 | #define HIGH(x) (x >> (BN_DIGIT_SIZE / 2)) | |
729 | #define L2H(x) (x << (BN_DIGIT_SIZE / 2)) | |
730 | ||
731 | ||
732 | BigNumDigit BnnDivideDigit (qq, nn, nl, d) | |
733 | ||
734 | register BigNum qq, nn; | |
735 | register int nl; | |
736 | BigNumDigit d; | |
737 | ||
738 | /* Performs the quotient: N div d => Q | |
739 | * Returns R = N mod d | |
740 | * Assumes leading digit of N < d, and d > 0. | |
741 | */ | |
742 | ||
743 | { | |
744 | /* test computed at compile time */ | |
745 | if (sizeof (BigNumProduct) > sizeof (BigNumDigit)) | |
746 | { | |
747 | register BigNumProduct quad; | |
748 | ||
749 | ||
750 | nn += nl; | |
751 | nl--; | |
752 | qq += nl; | |
753 | quad = *(--nn); | |
754 | ||
755 | while (nl != 0) | |
756 | { | |
757 | nl--; | |
758 | quad = (quad << BN_DIGIT_SIZE) | *(--nn); | |
759 | *(--qq) = quad / d; | |
760 | quad = quad % d; | |
761 | } | |
762 | ||
763 | return (quad); | |
764 | } | |
765 | else | |
766 | { | |
767 | int k; | |
768 | int orig_nl; | |
769 | BigNumDigit rh; /* Two halves of current remainder */ | |
770 | BigNumDigit rl; /* Correspond to quad above */ | |
771 | register BigNumDigit qa; /* Current appr. to quotient */ | |
772 | register BigNumDigit ph, pl; /* product of c and qa */ | |
773 | BigNumDigit ch, cl, prev_qq; | |
774 | ||
775 | ||
776 | /* Normalize divisor */ | |
777 | k = BnnNumLeadingZeroBitsInDigit (d); | |
778 | if (k != 0) | |
779 | { | |
780 | prev_qq = qq[-1]; | |
781 | orig_nl = nl; | |
782 | d <<= k; | |
783 | BnnShiftLeft (nn, nl, k); | |
784 | } | |
785 | ||
786 | nn += nl; | |
787 | nl--; | |
788 | qq += nl; | |
789 | ||
790 | ch = HIGH (d); | |
791 | cl = LOW (d); | |
792 | ||
793 | rl = *(--nn); | |
794 | ||
795 | while (nl != 0) | |
796 | { | |
797 | nl--; | |
798 | rh = rl; | |
799 | rl = *(--nn); | |
800 | qa = rh / ch; /* appr. quotient */ | |
801 | ||
802 | /* Compute ph, pl */ | |
803 | pl = cl * qa; | |
804 | ph = ch * qa; | |
805 | ph += HIGH (pl); | |
806 | pl = L2H (pl); | |
807 | ||
808 | /* While ph:pl > rh:rl, decrement qa, adjust qh:ql */ | |
809 | while (ph > rh || ph == rh && pl > rl) | |
810 | { | |
811 | qa--; | |
812 | SUB (ph, pl, ch, L2H (cl)); | |
813 | } | |
814 | ||
815 | SUB (rh, rl, ph, pl); | |
816 | ||
817 | /* Top half of quotient is correct; save it */ | |
818 | *(--qq) = L2H (qa); | |
819 | qa = (L2H (rh) | HIGH (rl)) / ch; | |
820 | ||
821 | /* Approx low half of q */ | |
822 | /* Compute ph, pl, again */ | |
823 | pl = cl * qa; | |
824 | ph = ch * qa; | |
825 | ph += HIGH (pl); | |
826 | pl = LOW (pl) | L2H (LOW (ph)); | |
827 | ph = HIGH (ph); | |
828 | ||
829 | /* While ph:pl > rh:rl, decrement qa, adjust qh:ql */ | |
830 | while (ph > rh || ph == rh && pl > rl) | |
831 | { | |
832 | qa--; | |
833 | SUB (ph, pl, 0, d); | |
834 | } | |
835 | ||
836 | /* Subtract ph:pl from rh:rl; we know rh will be 0 */ | |
837 | rl -= pl; | |
838 | *qq |= qa; | |
839 | } | |
840 | ||
841 | /* Denormalize dividend */ | |
842 | if (k != 0) { | |
843 | if((qq > nn) && (qq < &nn[orig_nl])) { | |
844 | /* Overlap between qq and nn. Care of *qq! */ | |
845 | orig_nl = (qq - nn); | |
846 | BnnShiftRight (nn, orig_nl, k); | |
847 | nn[orig_nl - 1] = prev_qq; | |
848 | } else if(qq == nn) { | |
849 | BnnShiftRight(&nn[orig_nl - 1], 1, k); | |
850 | } else { | |
851 | BnnShiftRight (nn, orig_nl, k); | |
852 | } } | |
853 | return (rl >> k); | |
854 | } | |
855 | } | |
856 | ||
857 | /***************************************/ | |
858 | ||
859 |