1 /* 128 bit integer arithmetic.
2  *
3  * Not optimized for speed.
4  *
5  * Copyright: Copyright D Language Foundation 2022.
6  * License:   $(LINK2 http://www.boost.org/LICENSE_1_0.txt, Boost License 1.0)
7  * Authors:   Walter Bright
8  * Source:    $(DRUNTIMESRC core/_int128.d)
9  */
10 
11 module core.int128;
12 
13 nothrow:
14 @safe:
15 @nogc:
16 
17 private alias I = long;
18 private alias U = ulong;
19 enum Ubits = uint(U.sizeof * 8);
20 
21 version (DigitalMars)
22 {
23     /* The alignment should follow target.stackAlign(),
24      * which is `isXmmSupported() ? 16 : (is64bit ? 8 : 4)
25      */
26     version (D_SIMD)
27         private enum Cent_alignment = 16;
28     else version (X86_64)
29         private enum Cent_alignment = 8;
30     else
31         private enum Cent_alignment = 4;
32 }
33 else
34 {
35     version (X86_64) private enum Cent_alignment = 16;
36     else             private enum Cent_alignment = (size_t.sizeof * 2);
37 }
38 
39 /**
40  * 128 bit integer type.
41  * See_also: $(REF Int128, std,int128).
42  */
43 align(Cent_alignment) struct Cent
44 {
45     version (LittleEndian)
46     {
47         U lo;  // low 64 bits
48         U hi;  // high 64 bits
49     }
50     else
51     {
52         U hi;  // high 64 bits
53         U lo;  // low 64 bits
54     }
55 }
56 
57 enum Cent One = { lo:1 };
58 enum Cent Zero = { lo:0 };
59 enum Cent MinusOne = neg(One);
60 
61 /*****************************
62  * Test against 0
63  * Params:
64  *      c = Cent to test
65  * Returns:
66  *      true if != 0
67  */
68 pure
69 bool tst(Cent c)
70 {
71     return c.hi || c.lo;
72 }
73 
74 
75 /*****************************
76  * Complement
77  * Params:
78  *      c = Cent to complement
79  * Returns:
80  *      complemented value
81  */
82 pure
83 Cent com(Cent c)
84 {
85     c.lo = ~c.lo;
86     c.hi = ~c.hi;
87     return c;
88 }
89 
90 /*****************************
91  * Negate
92  * Params:
93  *      c = Cent to negate
94  * Returns:
95  *      negated value
96  */
97 pure
98 Cent neg(Cent c)
99 {
100     if (c.lo == 0)
101         c.hi = -c.hi;
102     else
103     {
104         c.lo = -c.lo;
105         c.hi = ~c.hi;
106     }
107     return c;
108 }
109 
110 /*****************************
111  * Increment
112  * Params:
113  *      c = Cent to increment
114  * Returns:
115  *      incremented value
116  */
117 pure
118 Cent inc(Cent c)
119 {
120     return add(c, One);
121 }
122 
123 /*****************************
124  * Decrement
125  * Params:
126  *      c = Cent to decrement
127  * Returns:
128  *      incremented value
129  */
130 pure
131 Cent dec(Cent c)
132 {
133     return sub(c, One);
134 }
135 
136 /*****************************
137  * Shift left one bit
138  * Params:
139  *      c = Cent to shift
140  * Returns:
141  *      shifted value
142  */
143 pure
144 Cent shl1(Cent c)
145 {
146     c.hi = (c.hi << 1) | (cast(I)c.lo < 0);
147     c.lo <<= 1;
148     return c;
149 }
150 
151 /*****************************
152  * Unsigned shift right one bit
153  * Params:
154  *      c = Cent to shift
155  * Returns:
156  *      shifted value
157  */
158 pure
159 Cent shr1(Cent c)
160 {
161     c.lo = (c.lo >> 1) | ((c.hi & 1) << (Ubits - 1));
162     c.hi >>= 1;
163     return c;
164 }
165 
166 
167 /*****************************
168  * Arithmetic shift right one bit
169  * Params:
170  *      c = Cent to shift
171  * Returns:
172  *      shifted value
173  */
174 pure
175 Cent sar1(Cent c)
176 {
177     c.lo = (c.lo >> 1) | ((c.hi & 1) << (Ubits - 1));
178     c.hi = cast(I)c.hi >> 1;
179     return c;
180 }
181 
182 /*****************************
183  * Shift left n bits
184  * Params:
185  *      c = Cent to shift
186  *      n = number of bits to shift
187  * Returns:
188  *      shifted value
189  */
190 pure
191 Cent shl(Cent c, uint n)
192 {
193     if (n >= Ubits * 2)
194         return Zero;
195 
196     if (n >= Ubits)
197     {
198         c.hi = c.lo << (n - Ubits);
199         c.lo = 0;
200     }
201     else
202     {
203         c.hi = ((c.hi << n) | (c.lo >> (Ubits - n - 1) >> 1));
204         c.lo = c.lo << n;
205     }
206     return c;
207 }
208 
209 /*****************************
210  * Unsigned shift right n bits
211  * Params:
212  *      c = Cent to shift
213  *      n = number of bits to shift
214  * Returns:
215  *      shifted value
216  */
217 pure
218 Cent shr(Cent c, uint n)
219 {
220     if (n >= Ubits * 2)
221         return Zero;
222 
223     if (n >= Ubits)
224     {
225         c.lo = c.hi >> (n - Ubits);
226         c.hi = 0;
227     }
228     else
229     {
230         c.lo = ((c.lo >> n) | (c.hi << (Ubits - n - 1) << 1));
231         c.hi = c.hi >> n;
232     }
233     return c;
234 }
235 
236 /*****************************
237  * Arithmetic shift right n bits
238  * Params:
239  *      c = Cent to shift
240  *      n = number of bits to shift
241  * Returns:
242  *      shifted value
243  */
244 pure
245 Cent sar(Cent c, uint n)
246 {
247     const signmask = -(c.hi >> (Ubits - 1));
248     const signshift = (Ubits * 2) - n;
249     c = shr(c, n);
250 
251     // Sign extend all bits beyond the precision of Cent.
252     if (n >= Ubits * 2)
253     {
254         c.hi = signmask;
255         c.lo = signmask;
256     }
257     else if (signshift >= Ubits * 2)
258     {
259     }
260     else if (signshift >= Ubits)
261     {
262         c.hi &= ~(U.max << (signshift - Ubits));
263         c.hi |= signmask << (signshift - Ubits);
264     }
265     else
266     {
267         c.hi = signmask;
268         c.lo &= ~(U.max << signshift);
269         c.lo |= signmask << signshift;
270     }
271     return c;
272 }
273 
274 /*****************************
275  * Rotate left one bit
276  * Params:
277  *      c = Cent to rotate
278  * Returns:
279  *      rotated value
280  */
281 pure
282 Cent rol1(Cent c)
283 {
284     int carry = cast(I)c.hi < 0;
285 
286     c.hi = (c.hi << 1) | (cast(I)c.lo < 0);
287     c.lo = (c.lo << 1) | carry;
288     return c;
289 }
290 
291 /*****************************
292  * Rotate right one bit
293  * Params:
294  *      c = Cent to rotate
295  * Returns:
296  *      rotated value
297  */
298 pure
299 Cent ror1(Cent c)
300 {
301     int carry = c.lo & 1;
302     c.lo = (c.lo >> 1) | (cast(U)(c.hi & 1) << (Ubits - 1));
303     c.hi = (c.hi >> 1) | (cast(U)carry << (Ubits - 1));
304     return c;
305 }
306 
307 
308 /*****************************
309  * Rotate left n bits
310  * Params:
311  *      c = Cent to rotate
312  *      n = number of bits to rotate
313  * Returns:
314  *      rotated value
315  */
316 pure
317 Cent rol(Cent c, uint n)
318 {
319     n &= Ubits * 2 - 1;
320     Cent l = shl(c, n);
321     Cent r = shr(c, Ubits * 2 - n);
322     return or(l, r);
323 }
324 
325 /*****************************
326  * Rotate right n bits
327  * Params:
328  *      c = Cent to rotate
329  *      n = number of bits to rotate
330  * Returns:
331  *      rotated value
332  */
333 pure
334 Cent ror(Cent c, uint n)
335 {
336     n &= Ubits * 2 - 1;
337     Cent r = shr(c, n);
338     Cent l = shl(c, Ubits * 2 - n);
339     return or(r, l);
340 }
341 
342 /****************************
343  * And c1 & c2.
344  * Params:
345  *      c1 = operand 1
346  *      c2 = operand 2
347  * Returns:
348  *      c1 & c2
349  */
350 pure
351 Cent and(Cent c1, Cent c2)
352 {
353     const Cent ret = { lo:c1.lo & c2.lo, hi:c1.hi & c2.hi };
354     return ret;
355 }
356 
357 /****************************
358  * Or c1 | c2.
359  * Params:
360  *      c1 = operand 1
361  *      c2 = operand 2
362  * Returns:
363  *      c1 | c2
364  */
365 pure
366 Cent or(Cent c1, Cent c2)
367 {
368     const Cent ret = { lo:c1.lo | c2.lo, hi:c1.hi | c2.hi };
369     return ret;
370 }
371 
372 /****************************
373  * Xor c1 ^ c2.
374  * Params:
375  *      c1 = operand 1
376  *      c2 = operand 2
377  * Returns:
378  *      c1 ^ c2
379  */
380 pure
381 Cent xor(Cent c1, Cent c2)
382 {
383     const Cent ret = { lo:c1.lo ^ c2.lo, hi:c1.hi ^ c2.hi };
384     return ret;
385 }
386 
387 /****************************
388  * Add c1 to c2.
389  * Params:
390  *      c1 = operand 1
391  *      c2 = operand 2
392  * Returns:
393  *      c1 + c2
394  */
395 pure
396 Cent add(Cent c1, Cent c2)
397 {
398     U r = cast(U)(c1.lo + c2.lo);
399     const Cent ret = { lo:r, hi:cast(U)(c1.hi + c2.hi + (r < c1.lo)) };
400     return ret;
401 }
402 
403 /****************************
404  * Subtract c2 from c1.
405  * Params:
406  *      c1 = operand 1
407  *      c2 = operand 2
408  * Returns:
409  *      c1 - c2
410  */
411 pure
412 Cent sub(Cent c1, Cent c2)
413 {
414     return add(c1, neg(c2));
415 }
416 
417 /****************************
418  * Multiply c1 * c2.
419  * Params:
420  *      c1 = operand 1
421  *      c2 = operand 2
422  * Returns:
423  *      c1 * c2
424  */
425 pure
426 Cent mul(Cent c1, Cent c2)
427 {
428     enum mulmask = (1UL << (Ubits / 2)) - 1;
429     enum mulshift = Ubits / 2;
430 
431     // This algorithm splits the operands into 4 words, then computes and sums
432     // the partial products of each part.
433     const c2l0 = c2.lo & mulmask;
434     const c2l1 = c2.lo >> mulshift;
435     const c2h0 = c2.hi & mulmask;
436     const c2h1 = c2.hi >> mulshift;
437 
438     const c1l0 = c1.lo & mulmask;
439     U r0 = c1l0 * c2l0;
440     U r1 = c1l0 * c2l1 + (r0 >> mulshift);
441     U r2 = c1l0 * c2h0 + (r1 >> mulshift);
442     U r3 = c1l0 * c2h1 + (r2 >> mulshift);
443 
444     const c1l1 = c1.lo >> mulshift;
445     r1 = c1l1 * c2l0 + (r1 & mulmask);
446     r2 = c1l1 * c2l1 + (r2 & mulmask) + (r1 >> mulshift);
447     r3 = c1l1 * c2h0 + (r3 & mulmask) + (r2 >> mulshift);
448 
449     const c1h0 = c1.hi & mulmask;
450     r2 = c1h0 * c2l0 + (r2 & mulmask);
451     r3 = c1h0 * c2l1 + (r3 & mulmask) + (r2 >> mulshift);
452 
453     const c1h1 = c1.hi >> mulshift;
454     r3 = c1h1 * c2l0 + (r3 & mulmask);
455 
456     const Cent ret = { lo:(r0 & mulmask) + (r1 & mulmask) * (mulmask + 1),
457                        hi:(r2 & mulmask) + (r3 & mulmask) * (mulmask + 1) };
458     return ret;
459 }
460 
461 
462 /****************************
463  * Unsigned divide c1 / c2.
464  * Params:
465  *      c1 = dividend
466  *      c2 = divisor
467  * Returns:
468  *      quotient c1 / c2
469  */
470 pure
471 Cent udiv(Cent c1, Cent c2)
472 {
473     Cent modulus;
474     return udivmod(c1, c2, modulus);
475 }
476 
477 /****************************
478  * Unsigned divide c1 / c2. The remainder after division is stored to modulus.
479  * Params:
480  *      c1 = dividend
481  *      c2 = divisor
482  *      modulus = set to c1 % c2
483  * Returns:
484  *      quotient c1 / c2
485  */
486 pure
487 Cent udivmod(Cent c1, Cent c2, out Cent modulus)
488 {
489     //printf("udiv c1(%llx,%llx) c2(%llx,%llx)\n", c1.lo, c1.hi, c2.lo, c2.hi);
490     // Based on "Unsigned Doubleword Division" in Hacker's Delight
491     import core.bitop;
492 
493     // Divides a 128-bit dividend by a 64-bit divisor.
494     // The result must fit in 64 bits.
495     static U udivmod128_64(Cent c1, U c2, out U modulus)
496     {
497         // We work in base 2^^32
498         enum base = 1UL << 32;
499         enum divmask = (1UL << (Ubits / 2)) - 1;
500         enum divshift = Ubits / 2;
501 
502         // Check for overflow and divide by 0
503         if (c1.hi >= c2)
504         {
505             modulus = 0UL;
506             return ~0UL;
507         }
508 
509         // Computes [num1 num0] / den
510         static uint udiv96_64(U num1, uint num0, U den)
511         {
512             // Extract both digits of the denominator
513             const den1 = cast(uint)(den >> divshift);
514             const den0 = cast(uint)(den & divmask);
515             // Estimate ret as num1 / den1, and then correct it
516             U ret = num1 / den1;
517             const t2 = (num1 % den1) * base + num0;
518             const t1 = ret * den0;
519             if (t1 > t2)
520                 ret -= (t1 - t2 > den) ? 2 : 1;
521             return cast(uint)ret;
522         }
523 
524         // Determine the normalization factor. We multiply c2 by this, so that its leading
525         // digit is at least half base. In binary this means just shifting left by the number
526         // of leading zeros, so that there's a 1 in the MSB.
527         // We also shift number by the same amount. This cannot overflow because c1.hi < c2.
528         const shift = (Ubits - 1) - bsr(c2);
529         c2 <<= shift;
530         U num2 = c1.hi;
531         num2 <<= shift;
532         num2 |= (c1.lo >> (-shift & 63)) & (-cast(I)shift >> 63);
533         c1.lo <<= shift;
534 
535         // Extract the low digits of the numerator (after normalizing)
536         const num1 = cast(uint)(c1.lo >> divshift);
537         const num0 = cast(uint)(c1.lo & divmask);
538 
539         // Compute q1 = [num2 num1] / c2
540         const q1 = udiv96_64(num2, num1, c2);
541         // Compute the true (partial) remainder
542         const rem = num2 * base + num1 - q1 * c2;
543         // Compute q0 = [rem num0] / c2
544         const q0 = udiv96_64(rem, num0, c2);
545 
546         modulus = (rem * base + num0 - q0 * c2) >> shift;
547         return (cast(U)q1 << divshift) | q0;
548     }
549 
550     // Special cases
551     if (!tst(c2))
552     {
553         // Divide by zero
554         modulus = Zero;
555         return com(modulus);
556     }
557     if (c1.hi == 0 && c2.hi == 0)
558     {
559         // Single precision divide
560         const Cent rem = { lo:c1.lo % c2.lo };
561         modulus = rem;
562         const Cent ret = { lo:c1.lo / c2.lo };
563         return ret;
564     }
565     if (c1.hi == 0)
566     {
567         // Numerator is smaller than the divisor
568         modulus = c1;
569         return Zero;
570     }
571     if (c2.hi == 0)
572     {
573         // Divisor is a 64-bit value, so we just need one 128/64 division.
574         // If c1 / c2 would overflow, break c1 up into two halves.
575         const q1 = (c1.hi < c2.lo) ? 0 : (c1.hi / c2.lo);
576         if (q1)
577             c1.hi = c1.hi % c2.lo;
578         Cent rem;
579         const q0 = udivmod128_64(c1, c2.lo, rem.lo);
580         modulus = rem;
581         const Cent ret = { lo:q0, hi:q1 };
582         return ret;
583     }
584 
585     // Full cent precision division.
586     // Here c2 >= 2^^64
587     // We know that c2.hi != 0, so count leading zeros is OK
588     // We have 0 <= shift <= 63
589     const shift = (Ubits - 1) - bsr(c2.hi);
590 
591     // Normalize the divisor so its MSB is 1
592     // v1 = (c2 << shift) >> 64
593     U v1 = shl(c2, shift).hi;
594 
595     // To ensure no overflow.
596     Cent u1 = shr1(c1);
597 
598     // Get quotient from divide unsigned operation.
599     U rem_ignored;
600     const Cent q1 = { lo:udivmod128_64(u1, v1, rem_ignored) };
601 
602     // Undo normalization and division of c1 by 2.
603     Cent quotient = shr(shl(q1, shift), 63);
604 
605     // Make quotient correct or too small by 1
606     if (tst(quotient))
607         quotient = dec(quotient);
608 
609     // Now quotient is correct.
610     // Compute rem = c1 - (quotient * c2);
611     Cent rem = sub(c1, mul(quotient, c2));
612 
613     // Check if remainder is larger than the divisor
614     if (uge(rem, c2))
615     {
616         // Increment quotient
617         quotient = inc(quotient);
618         // Subtract c2 from remainder
619         rem = sub(rem, c2);
620     }
621     modulus = rem;
622     //printf("quotient "); print(quotient);
623     //printf("modulus  "); print(modulus);
624     return quotient;
625 }
626 
627 
628 /****************************
629  * Signed divide c1 / c2.
630  * Params:
631  *      c1 = dividend
632  *      c2 = divisor
633  * Returns:
634  *      quotient c1 / c2
635  */
636 pure
637 Cent div(Cent c1, Cent c2)
638 {
639     Cent modulus;
640     return divmod(c1, c2, modulus);
641 }
642 
643 /****************************
644  * Signed divide c1 / c2. The remainder after division is stored to modulus.
645  * Params:
646  *      c1 = dividend
647  *      c2 = divisor
648  *      modulus = set to c1 % c2
649  * Returns:
650  *      quotient c1 / c2
651  */
652 pure
653 Cent divmod(Cent c1, Cent c2, out Cent modulus)
654 {
655     /* Muck about with the signs so we can use the unsigned divide
656      */
657     if (cast(I)c1.hi < 0)
658     {
659         if (cast(I)c2.hi < 0)
660         {
661             Cent r = udivmod(neg(c1), neg(c2), modulus);
662             modulus = neg(modulus);
663             return r;
664         }
665         Cent r = neg(udivmod(neg(c1), c2, modulus));
666         modulus = neg(modulus);
667         return r;
668     }
669     else if (cast(I)c2.hi < 0)
670     {
671         return neg(udivmod(c1, neg(c2), modulus));
672     }
673     else
674         return udivmod(c1, c2, modulus);
675 }
676 
677 /****************************
678  * If c1 > c2 unsigned
679  * Params:
680  *      c1 = operand 1
681  *      c2 = operand 2
682  * Returns:
683  *      true if c1 > c2
684  */
685 pure
686 bool ugt(Cent c1, Cent c2)
687 {
688     return (c1.hi == c2.hi) ? (c1.lo > c2.lo) : (c1.hi > c2.hi);
689 }
690 
691 /****************************
692  * If c1 >= c2 unsigned
693  * Params:
694  *      c1 = operand 1
695  *      c2 = operand 2
696  * Returns:
697  *      true if c1 >= c2
698  */
699 pure
700 bool uge(Cent c1, Cent c2)
701 {
702     return !ugt(c2, c1);
703 }
704 
705 /****************************
706  * If c1 < c2 unsigned
707  * Params:
708  *      c1 = operand 1
709  *      c2 = operand 2
710  * Returns:
711  *      true if c1 < c2
712  */
713 pure
714 bool ult(Cent c1, Cent c2)
715 {
716     return ugt(c2, c1);
717 }
718 
719 /****************************
720  * If c1 <= c2 unsigned
721  * Params:
722  *      c1 = operand 1
723  *      c2 = operand 2
724  * Returns:
725  *      true if c1 <= c2
726  */
727 pure
728 bool ule(Cent c1, Cent c2)
729 {
730     return !ugt(c1, c2);
731 }
732 
733 /****************************
734  * If c1 > c2 signed
735  * Params:
736  *      c1 = operand 1
737  *      c2 = operand 2
738  * Returns:
739  *      true if c1 > c2
740  */
741 pure
742 bool gt(Cent c1, Cent c2)
743 {
744     return (c1.hi == c2.hi)
745         ? (c1.lo > c2.lo)
746         : (cast(I)c1.hi > cast(I)c2.hi);
747 }
748 
749 /****************************
750  * If c1 >= c2 signed
751  * Params:
752  *      c1 = operand 1
753  *      c2 = operand 2
754  * Returns:
755  *      true if c1 >= c2
756  */
757 pure
758 bool ge(Cent c1, Cent c2)
759 {
760     return !gt(c2, c1);
761 }
762 
763 /****************************
764  * If c1 < c2 signed
765  * Params:
766  *      c1 = operand 1
767  *      c2 = operand 2
768  * Returns:
769  *      true if c1 < c2
770  */
771 pure
772 bool lt(Cent c1, Cent c2)
773 {
774     return gt(c2, c1);
775 }
776 
777 /****************************
778  * If c1 <= c2 signed
779  * Params:
780  *      c1 = operand 1
781  *      c2 = operand 2
782  * Returns:
783  *      true if c1 <= c2
784  */
785 pure
786 bool le(Cent c1, Cent c2)
787 {
788     return !gt(c1, c2);
789 }
790 
791 /*******************************************************/
792 
793 version (unittest)
794 {
795     version (none)
796     {
797         import core.stdc.stdio;
798 
799         @trusted
800         void print(Cent c)
801         {
802             printf("%lld, %lld\n", cast(ulong)c.lo, cast(ulong)c.hi);
803             printf("x%llx, x%llx\n", cast(ulong)c.lo, cast(ulong)c.hi);
804         }
805     }
806 }
807 
808 unittest
809 {
810     const Cent C0 = Zero;
811     const Cent C1 = One;
812     const Cent C2 = { lo:2 };
813     const Cent C3 = { lo:3 };
814     const Cent C5 = { lo:5 };
815     const Cent C10 = { lo:10 };
816     const Cent C20 = { lo:20 };
817     const Cent C30 = { lo:30 };
818     const Cent C100 = { lo:100 };
819 
820     const Cent Cm1 =  neg(One);
821     const Cent Cm3 =  neg(C3);
822     const Cent Cm10 = neg(C10);
823 
824     const Cent C3_1 = { lo:1, hi:3 };
825     const Cent C3_2 = { lo:2, hi:3 };
826     const Cent C4_8  = { lo:8, hi:4 };
827     const Cent C5_0  = { lo:0, hi:5 };
828     const Cent C7_1 = { lo:1, hi:7 };
829     const Cent C7_9 = { lo:9, hi:7 };
830     const Cent C9_3 = { lo:3, hi:9 };
831     const Cent C10_0 = { lo:0, hi:10 };
832     const Cent C10_1 = { lo:1, hi:10 };
833     const Cent C10_3 = { lo:3, hi:10 };
834     const Cent C11_3 = { lo:3, hi:11 };
835     const Cent C20_0 = { lo:0, hi:20 };
836     const Cent C90_30 = { lo:30, hi:90 };
837 
838     const Cent Cm10_0 = inc(com(C10_0)); // Cent(lo=0,  hi=-10);
839     const Cent Cm10_1 = inc(com(C10_1)); // Cent(lo=-1, hi=-11);
840     const Cent Cm10_3 = inc(com(C10_3)); // Cent(lo=-3, hi=-11);
841     const Cent Cm20_0 = inc(com(C20_0)); // Cent(lo=0,  hi=-20);
842 
843     enum Cent Cs_3 = { lo:3, hi:I.min };
844 
845     const Cent Cbig_1 = { lo:0xa3ccac1832952398, hi:0xc3ac542864f652f8 };
846     const Cent Cbig_2 = { lo:0x5267b85f8a42fc20, hi:0 };
847     const Cent Cbig_3 = { lo:0xf0000000ffffffff, hi:0 };
848 
849     /************************/
850 
851     assert( ugt(C1, C0) );
852     assert( ult(C1, C2) );
853     assert( uge(C1, C0) );
854     assert( ule(C1, C2) );
855 
856     assert( !ugt(C0, C1) );
857     assert( !ult(C2, C1) );
858     assert( !uge(C0, C1) );
859     assert( !ule(C2, C1) );
860 
861     assert( !ugt(C1, C1) );
862     assert( !ult(C1, C1) );
863     assert( uge(C1, C1) );
864     assert( ule(C2, C2) );
865 
866     assert( ugt(C10_3, C10_1) );
867     assert( ugt(C11_3, C10_3) );
868     assert( !ugt(C9_3, C10_3) );
869     assert( !ugt(C9_3, C9_3) );
870 
871     assert( gt(C2, C1) );
872     assert( !gt(C1, C2) );
873     assert( !gt(C1, C1) );
874     assert( gt(C0, Cm1) );
875     assert( gt(Cm1, neg(C10)));
876     assert( !gt(Cm1, Cm1) );
877     assert( !gt(Cm1, C0) );
878 
879     assert( !lt(C2, C1) );
880     assert( !le(C2, C1) );
881     assert( ge(C2, C1) );
882 
883     assert(neg(C10_0) == Cm10_0);
884     assert(neg(C10_1) == Cm10_1);
885     assert(neg(C10_3) == Cm10_3);
886 
887     assert(add(C7_1,C3_2) == C10_3);
888     assert(sub(C1,C2) == Cm1);
889 
890     assert(inc(C3_1) == C3_2);
891     assert(dec(C3_2) == C3_1);
892 
893     assert(shl(C10,0) == C10);
894     assert(shl(C10,Ubits) == C10_0);
895     assert(shl(C10,1) == C20);
896     assert(shl(C10,Ubits * 2) == C0);
897     assert(shr(C10_0,0) == C10_0);
898     assert(shr(C10_0,Ubits) == C10);
899     assert(shr(C10_0,Ubits - 1) == C20);
900     assert(shr(C10_0,Ubits + 1) == C5);
901     assert(shr(C10_0,Ubits * 2) == C0);
902     assert(sar(C10_0,0) == C10_0);
903     assert(sar(C10_0,Ubits) == C10);
904     assert(sar(C10_0,Ubits - 1) == C20);
905     assert(sar(C10_0,Ubits + 1) == C5);
906     assert(sar(C10_0,Ubits * 2) == C0);
907     assert(sar(Cm1,Ubits * 2) == Cm1);
908 
909     assert(shl1(C10) == C20);
910     assert(shr1(C10_0) == C5_0);
911     assert(sar1(C10_0) == C5_0);
912     assert(sar1(Cm1) == Cm1);
913 
914     Cent modulus;
915 
916     assert(udiv(C10,C2) == C5);
917     assert(udivmod(C10,C2, modulus) ==  C5);   assert(modulus == C0);
918     assert(udivmod(C10,C3, modulus) ==  C3);   assert(modulus == C1);
919     assert(udivmod(C10,C0, modulus) == Cm1);   assert(modulus == C0);
920     assert(udivmod(C2,C90_30, modulus) == C0); assert(modulus == C2);
921     assert(udiv(mul(C90_30, C2), C2) == C90_30);
922     assert(udiv(mul(C90_30, C2), C90_30) == C2);
923 
924     assert(div(C10,C3) == C3);
925     assert(divmod( C10,  C3, modulus) ==  C3); assert(modulus ==  C1);
926     assert(divmod(Cm10,  C3, modulus) == Cm3); assert(modulus == Cm1);
927     assert(divmod( C10, Cm3, modulus) == Cm3); assert(modulus ==  C1);
928     assert(divmod(Cm10, Cm3, modulus) ==  C3); assert(modulus == Cm1);
929     assert(divmod(C2, C90_30, modulus) == C0); assert(modulus == C2);
930     assert(div(mul(C90_30, C2), C2) == C90_30);
931     assert(div(mul(C90_30, C2), C90_30) == C2);
932 
933     const Cent Cb1divb2 = { lo:0x4496aa309d4d4a2f, hi:U.max };
934     const Cent Cb1modb2 = { lo:0xd83203d0fdc799b8, hi:U.max };
935     assert(divmod(Cbig_1, Cbig_2, modulus) == Cb1divb2);
936     assert(modulus == Cb1modb2);
937 
938     const Cent Cb1udivb2 = { lo:0x5fe0e9bace2bedad, hi:2 };
939     const Cent Cb1umodb2 = { lo:0x2c923125a68721f8, hi:0 };
940     assert(udivmod(Cbig_1, Cbig_2, modulus) == Cb1udivb2);
941     assert(modulus == Cb1umodb2);
942 
943     const Cent Cb1divb3 = { lo:0xbfa6c02b5aff8b86, hi:U.max };
944     const Cent Cb1udivb3 = { lo:0xd0b7d13b48cb350f, hi:0 };
945     assert(div(Cbig_1, Cbig_3) == Cb1divb3);
946     assert(udiv(Cbig_1, Cbig_3) == Cb1udivb3);
947 
948     assert(mul(Cm10, C1) == Cm10);
949     assert(mul(C1, Cm10) == Cm10);
950     assert(mul(C9_3, C10) == C90_30);
951     assert(mul(Cs_3, C10) == C30);
952     assert(mul(Cm10, Cm10) == C100);
953     assert(mul(C20_0, Cm1) == Cm20_0);
954 
955     assert( or(C4_8, C3_1) == C7_9);
956     assert(and(C4_8, C7_9) == C4_8);
957     assert(xor(C4_8, C7_9) == C3_1);
958 
959     assert(rol(Cm1,  1) == Cm1);
960     assert(ror(Cm1, 45) == Cm1);
961     assert(rol(ror(C7_9, 5), 5) == C7_9);
962     assert(rol(C7_9, 1) == rol1(C7_9));
963     assert(ror(C7_9, 1) == ror1(C7_9));
964 }