*** /dev/null Sat Nov 23 01:27:14 2024
--- - Sat Nov 23 01:27:29 2024
***************
*** 0 ****
--- 1,583 ----
+ /*
+ ** OSSP ui64 - 64-Bit Arithmetic
+ ** Copyright (c) 2002 Ralf S. Engelschall <rse@engelschall.com>
+ ** Copyright (c) 2002 The OSSP Project <http://www.ossp.org/>
+ ** Copyright (c) 2002 Cable & Wireless Deutschland <http://www.cw.com/de/>
+ **
+ ** This file is part of OSSP ui64, a 64-bit arithmetic library
+ ** which can be found at http://www.ossp.org/pkg/lib/ui64/.
+ **
+ ** Permission to use, copy, modify, and distribute this software for
+ ** any purpose with or without fee is hereby granted, provided that
+ ** the above copyright notice and this permission notice appear in all
+ ** copies.
+ **
+ ** THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESSED OR IMPLIED
+ ** WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
+ ** MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
+ ** IN NO EVENT SHALL THE AUTHORS AND COPYRIGHT HOLDERS AND THEIR
+ ** CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+ ** SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+ ** LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
+ ** USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
+ ** ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+ ** OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
+ ** OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ ** SUCH DAMAGE.
+ **
+ ** ui64.c: implementation of 64-bit unsigned integer arithmetic
+ */
+
+ #include <string.h>
+ #include <ctype.h>
+
+ #include "tai_ui64.h"
+
+ #define UI64_BASE 256 /* 2^8 */
+ #define UI64_DIGITS 8 /* 8*8 = 64 bit */
+ #define UIXX_T(n) struct { unsigned char x[n]; }
+
+ /* fill an ui64_t with a sequence of a particular digit */
+ #define ui64_fill(__x, __n) \
+ do { int __i; \
+ for (__i = 0; __i < UI64_DIGITS; __i++) \
+ (__x).x[__i] = (__n); \
+ } while(0)
+
+ /* the value zero */
+ ui64_t ui64_zero(void)
+ {
+ ui64_t z;
+
+ ui64_fill(z, 0);
+ return z;
+ }
+
+ /* the maximum value */
+ ui64_t ui64_max(void)
+ {
+ ui64_t z;
+
+ ui64_fill(z, UI64_BASE-1);
+ return z;
+ }
+
+ /* convert ISO-C "unsigned long" into internal format */
+ ui64_t ui64_n2i(unsigned long n)
+ {
+ ui64_t z;
+ int i;
+
+ i = 0;
+ do {
+ z.x[i++] = (n % UI64_BASE);
+ } while ((n /= UI64_BASE) > 0 && i < UI64_DIGITS);
+ for ( ; i < UI64_DIGITS; i++)
+ z.x[i] = 0;
+ return z;
+ }
+
+ /* convert internal format into ISO-C "unsigned long";
+ truncates if sizeof(unsigned long) is less than UI64_DIGITS! */
+ unsigned long ui64_i2n(ui64_t x)
+ {
+ unsigned long n;
+ int i;
+
+ n = 0;
+ i = (int)sizeof(n);
+ if (i > UI64_DIGITS)
+ i = UI64_DIGITS;
+ while (--i >= 0) {
+ n = (n * UI64_BASE) + x.x[i];
+ }
+ return n;
+ }
+
+ /* convert string representation of arbitrary base into internal format */
+ ui64_t ui64_s2i(const char *str, char **end, int base)
+ {
+ ui64_t z;
+ const char *cp;
+ int carry;
+ static char map[] = {
+ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, /* 0...9 */
+ 36, 36, 36, 36, 36, 36, 36, /* illegal chars */
+ 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, /* A...M */
+ 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, /* N...Z */
+ 36, 36, 36, 36, 36, 36, /* illegal chars */
+ 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, /* a...m */
+ 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35 /* m...z */
+ };
+
+ ui64_fill(z, 0);
+ if (str == NULL || (base < 2 || base > 36))
+ return z;
+ cp = str;
+ while (*cp != '\0' && isspace((int)(*cp)))
+ cp++;
+ while ( *cp != '\0'
+ && isalnum((int)(*cp))
+ && map[(int)(*cp)-'0'] < base) {
+ z = ui64_muln(z, base, &carry);
+ if (carry)
+ break;
+ z = ui64_addn(z, map[(int)(*cp)-'0'], &carry);
+ if (carry)
+ break;
+ cp++;
+ }
+ if (end != NULL)
+ *end = (char *)cp;
+ return z;
+ }
+
+ /* convert internal format into string representation of arbitrary base */
+ char *ui64_i2s(ui64_t x, char *str, size_t len, int base)
+ {
+ static char map[] = "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZ";
+ char c;
+ int r;
+ int n;
+ int i, j;
+
+ if (str == NULL || len < 2 || (base < 2 || base > 36))
+ return NULL;
+ n = ui64_len(x);
+ i = 0;
+ do {
+ x = ui64_divn(x, base, &r);
+ str[i++] = map[r];
+ while (n > 1 && x.x[n-1] == 0)
+ n--;
+ } while (i < (len-1) && (n > 1 || x.x[0] != 0));
+ str[i] = '\0';
+ for (j = 0; j < --i; j++) {
+ c = str[j];
+ str[j] = str[i];
+ str[i] = c;
+ }
+ return str;
+ }
+
+ /* addition of two ui64_t */
+ ui64_t ui64_add(ui64_t x, ui64_t y, ui64_t *ov)
+ {
+ ui64_t z;
+ int carry;
+ int i;
+
+ carry = 0;
+ for (i = 0; i < UI64_DIGITS; i++) {
+ carry += (x.x[i] + y.x[i]);
+ z.x[i] = (carry % UI64_BASE);
+ carry /= UI64_BASE;
+ }
+ if (ov != NULL)
+ *ov = ui64_n2i((unsigned long)carry);
+ return z;
+ }
+
+ /* addition of an ui64_t and a single digit */
+ ui64_t ui64_addn(ui64_t x, int y, int *ov)
+ {
+ ui64_t z;
+ int i;
+
+ for (i = 0; i < UI64_DIGITS; i++) {
+ y += x.x[i];
+ z.x[i] = (y % UI64_BASE);
+ y /= UI64_BASE;
+ }
+ if (ov != NULL)
+ *ov = y;
+ return z;
+ }
+
+ /* subtraction of two ui64_t */
+ ui64_t ui64_sub(ui64_t x, ui64_t y, ui64_t *ov)
+ {
+ ui64_t z;
+ int borrow;
+ int i;
+ int d;
+
+ borrow = 0;
+ for (i = 0; i < UI64_DIGITS; i++) {
+ d = ((x.x[i] + UI64_BASE) - borrow - y.x[i]);
+ z.x[i] = (d % UI64_BASE);
+ borrow = (1 - (d/UI64_BASE));
+ }
+ if (ov != NULL)
+ *ov = ui64_n2i((unsigned long)borrow);
+ return z;
+ }
+
+ /* subtraction of an ui64_t and a single digit */
+ ui64_t ui64_subn(ui64_t x, int y, int *ov)
+ {
+ ui64_t z;
+ int i;
+ int d;
+
+ for (i = 0; i < UI64_DIGITS; i++) {
+ d = (x.x[i] + UI64_BASE) - y;
+ z.x[i] = (d % UI64_BASE);
+ y = (1 - (d/UI64_BASE));
+ }
+ if (ov != NULL)
+ *ov = y;
+ return z;
+ }
+
+ /*
+ 7 3 2
+ * 9 4 2 8
+ ---------
+ 5 8 5 6
+ + 1 4 6 4
+ + 2 9 2 8
+ + 6 5 8 8
+ ---------------
+ = 6 9 0 1 2 9 6
+ */
+
+ ui64_t ui64_mul(ui64_t x, ui64_t y, ui64_t *ov)
+ {
+ UIXX_T(UI64_DIGITS+UI64_DIGITS) zx;
+ ui64_t z;
+ int carry;
+ int i, j;
+
+ /* clear temporary result buffer */
+ for (i = 0; i < (UI64_DIGITS+UI64_DIGITS); i++)
+ zx.x[i] = 0;
+
+ /* perform multiplication operation */
+ for (i = 0; i < UI64_DIGITS; i++) {
+ /* calculate partial product and immediately add to z */
+ carry = 0;
+ for (j = 0; j < UI64_DIGITS; j++) {
+ carry += (x.x[i] * y.x[j]) + zx.x[i+j];
+ zx.x[i+j] = (carry % UI64_BASE);
+ carry /= UI64_BASE;
+ }
+ /* add carry to remaining digits in z */
+ for ( ; j < UI64_DIGITS + UI64_DIGITS - i; j++) {
+ carry += zx.x[i+j];
+ zx.x[i+j] = (carry % UI64_BASE);
+ carry /= UI64_BASE;
+ }
+ }
+
+ /* provide result by splitting zx into z and ov */
+ memcpy(z.x, zx.x, UI64_DIGITS);
+ if (ov != NULL)
+ memcpy((*ov).x, &zx.x[UI64_DIGITS], UI64_DIGITS);
+
+ return z;
+ }
+
+ ui64_t ui64_muln(ui64_t x, int y, int *ov)
+ {
+ ui64_t z;
+ int carry;
+ int i;
+
+ carry = 0;
+ for (i = 0; i < UI64_DIGITS; i++) {
+ carry += (x.x[i] * y);
+ z.x[i] = (carry % UI64_BASE);
+ carry /= UI64_BASE;
+ }
+ if (ov != NULL)
+ *ov = carry;
+ return z;
+ }
+
+ /*
+ = 2078 [q]
+ 0615367 [x] : 296 [y]
+ -0592 [dq]
+ -----
+ = 0233
+ -0000 [dq]
+ -----
+ = 2336
+ -2072 [dq]
+ -----
+ = 2647
+ -2308 [dq]
+ -----
+ = 279 [r]
+ */
+ ui64_t ui64_div(ui64_t x, ui64_t y, ui64_t *ov)
+ {
+ ui64_t q;
+ ui64_t r;
+ int i;
+ int n, m;
+ int ovn;
+
+ /* determine actual number of involved digits */
+ n = ui64_len(x);
+ m = ui64_len(y);
+
+ if (m == 1) {
+ /* simple case #1: reduceable to ui64_divn() */
+ if (y.x[0] == 0) {
+ /* error case: division by zero! */
+ ui64_fill(q, 0);
+ ui64_fill(r, 0);
+ }
+ else {
+ q = ui64_divn(x, y.x[0], &ovn);
+ ui64_fill(r, 0);
+ r.x[0] = (unsigned char)ovn;
+ }
+
+ } else if (n < m) {
+ /* simple case #2: everything is in the remainder */
+ ui64_fill(q, 0);
+ r = x;
+
+ } else { /* n >= m, m > 1 */
+ /* standard case: x[0..n] / y[0..m] */
+ UIXX_T(UI64_DIGITS+1) rx;
+ UIXX_T(UI64_DIGITS+1) dq;
+ ui64_t t;
+ int km;
+ int k;
+ int qk;
+ unsigned long y2;
+ unsigned long r3;
+ int borrow;
+ int d;
+
+ /* rx is x with a leading zero in order to make
+ sure that n > m and not just n >= m */
+ memcpy(rx.x, x.x, UI64_DIGITS);
+ rx.x[UI64_DIGITS] = 0;
+
+ for (k = n - m; k >= 0; k--) {
+ /* efficiently compute qk by guessing
+ qk := rx[k+m-2...k+m]/y[m-2...m-1] */
+ km = k + m;
+ y2 = (y.x[m-1]*UI64_BASE) + y.x[m-2];
+ r3 = (rx.x[km]*(UI64_BASE*UI64_BASE)) +
+ (rx.x[km-1]*UI64_BASE) + rx.x[km-2];
+ qk = r3 / y2;
+ if (qk >= UI64_BASE)
+ qk = UI64_BASE - 1;
+
+ /* dq := y*qk (post-adjust qk if guessed incorrectly) */
+ t = ui64_muln(y, qk, &ovn);
+ memcpy(dq.x, t.x, UI64_DIGITS);
+ dq.x[m] = (unsigned char)ovn;
+ for (i = m; i > 0; i--)
+ if (rx.x[i+k] != dq.x[i])
+ break;
+ if (rx.x[i+k] < dq.x[i]) {
+ t = ui64_muln(y, --qk, &ovn);
+ memcpy(dq.x, t.x, UI64_DIGITS);
+ dq.x[m] = (unsigned char)ovn;
+ }
+
+ /* store qk */
+ q.x[k] = (unsigned char)qk;
+
+ /* rx := rx - dq*(b^k) */
+ borrow = 0;
+ for (i = 0; i < m+1; i++) {
+ d = ((rx.x[k+i] + UI64_BASE) - borrow - dq.x[i]);
+ rx.x[k+i] = (d % UI64_BASE);
+ borrow = (1 - (d/UI64_BASE));
+ }
+ }
+ memcpy(r.x, rx.x, m);
+
+ /* fill out results with leading zeros */
+ for (i = n-m+1; i < UI64_DIGITS; i++)
+ q.x[i] = 0;
+ for (i = m; i < UI64_DIGITS; i++)
+ r.x[i] = 0;
+ }
+
+ /* provide results */
+ if (ov != NULL)
+ *ov = r;
+ return q;
+ }
+
+ ui64_t ui64_divn(ui64_t x, int y, int *ov)
+ {
+ ui64_t z;
+ unsigned int carry;
+ int i;
+
+ carry = 0;
+ for (i = (UI64_DIGITS - 1); i >= 0; i--) {
+ carry = (carry * UI64_BASE) + x.x[i];
+ z.x[i] = (carry / y);
+ carry %= y;
+ }
+ if (ov != NULL)
+ *ov = carry;
+ return z;
+ }
+
+ ui64_t ui64_and(ui64_t x, ui64_t y)
+ {
+ ui64_t z;
+ int i;
+
+ for (i = 0; i < UI64_DIGITS; i++)
+ z.x[i] = (x.x[i] & y.x[i]);
+ return z;
+ }
+
+ ui64_t ui64_or(ui64_t x, ui64_t y)
+ {
+ ui64_t z;
+ int i;
+
+ for (i = 0; i < UI64_DIGITS; i++)
+ z.x[i] = (x.x[i] | y.x[i]);
+ return z;
+ }
+
+ ui64_t ui64_xor(ui64_t x, ui64_t y)
+ {
+ ui64_t z;
+ int i;
+
+ for (i = 0; i < UI64_DIGITS; i++)
+ z.x[i] = ((x.x[i] & ~(y.x[i])) | (~(x.x[i]) & (y.x[i])));
+ return z;
+ }
+
+ ui64_t ui64_not(ui64_t x)
+ {
+ ui64_t z;
+ int i;
+
+ for (i = 0; i < UI64_DIGITS; i++)
+ z.x[i] = ~(x.x[i]);
+ return z;
+ }
+
+ ui64_t ui64_rol(ui64_t x, int s, ui64_t *ov)
+ {
+ UIXX_T(UI64_DIGITS+UI64_DIGITS) zx;
+ ui64_t z;
+ int i;
+ int carry;
+
+ if (s <= 0) {
+ /* no shift at all */
+ if (ov != NULL)
+ *ov = ui64_zero();
+ return x;
+ }
+ else if (s > 64) {
+ /* too large shift */
+ if (ov != NULL)
+ *ov = ui64_zero();
+ return ui64_zero();
+ }
+ else if (s == 64) {
+ /* maximum shift */
+ if (ov != NULL)
+ *ov = x;
+ return ui64_zero();
+ }
+ else { /* regular shift */
+ /* shift (logically) left by s/8 bytes */
+ for (i = 0; i < UI64_DIGITS+UI64_DIGITS; i++)
+ zx.x[i] = 0;
+ for (i = 0; i < UI64_DIGITS; i++)
+ zx.x[i+(s/8)] = x.x[i];
+ /* shift (logically) left by remaining s%8 bits */
+ s %= 8;
+ if (s > 0) {
+ carry = 0;
+ for (i = 0; i < UI64_DIGITS+UI64_DIGITS; i++) {
+ carry += (zx.x[i] * (1 << s));
+ zx.x[i] = (carry % UI64_BASE);
+ carry /= UI64_BASE;
+ }
+ }
+ memcpy(z.x, zx.x, UI64_DIGITS);
+ if (ov != NULL)
+ memcpy((*ov).x, &zx.x[UI64_DIGITS], UI64_DIGITS);
+ }
+ return z;
+ }
+
+ ui64_t ui64_ror(ui64_t x, int s, ui64_t *ov)
+ {
+ UIXX_T(UI64_DIGITS+UI64_DIGITS) zx;
+ ui64_t z;
+ int i;
+ int carry;
+
+ if (s <= 0) {
+ /* no shift at all */
+ if (ov != NULL)
+ *ov = ui64_zero();
+ return x;
+ }
+ else if (s > 64) {
+ /* too large shift */
+ if (ov != NULL)
+ *ov = ui64_zero();
+ return ui64_zero();
+ }
+ else if (s == 64) {
+ /* maximum shift */
+ if (ov != NULL)
+ *ov = x;
+ return ui64_zero();
+ }
+ else { /* regular shift */
+ /* shift (logically) right by s/8 bytes */
+ for (i = 0; i < UI64_DIGITS+UI64_DIGITS; i++)
+ zx.x[i] = 0;
+ for (i = 0; i < UI64_DIGITS; i++)
+ zx.x[UI64_DIGITS+i-(s/8)] = x.x[i];
+ /* shift (logically) right by remaining s%8 bits */
+ s %= 8;
+ if (s > 0) {
+ carry = 0;
+ for (i = (UI64_DIGITS+UI64_DIGITS - 1); i >= 0; i--) {
+ carry = (carry * UI64_BASE) + zx.x[i];
+ zx.x[i] = (carry / (1 << s));
+ carry %= (1 << s);
+ }
+ }
+ memcpy(z.x, &zx.x[UI64_DIGITS], UI64_DIGITS);
+ if (ov != NULL)
+ memcpy((*ov).x, zx.x, UI64_DIGITS);
+ }
+ return z;
+ }
+
+ int ui64_cmp(ui64_t x, ui64_t y)
+ {
+ int i;
+
+ i = UI64_DIGITS - 1;
+ while (i > 0 && x.x[i] == y.x[i])
+ i--;
+ return (x.x[i] - y.x[i]);
+ }
+
+ int ui64_len(ui64_t x)
+ {
+ int i;
+
+ for (i = UI64_DIGITS; i > 1 && x.x[i-1] == 0; i--)
+ ;
+ return i;
+ }
+
|
