summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
Diffstat
-rw-r--r--libtommath/bn_mp_div_3.c63
-rw-r--r--libtommath/bn_mp_expt_u32.c46
-rw-r--r--libtommath/bn_mp_log_u32.c180
-rw-r--r--libtommath/bn_mp_root_u32.c139
4 files changed, 0 insertions, 428 deletions
diff --git a/libtommath/bn_mp_div_3.c b/libtommath/bn_mp_div_3.c
deleted file mode 100644
index 3a23fdf..0000000
--- a/libtommath/bn_mp_div_3.c
+++ /dev/null
@@ -1,63 +0,0 @@
-#include "tommath_private.h"
-#ifdef BN_MP_DIV_3_C
-/* LibTomMath, multiple-precision integer library -- Tom St Denis */
-/* SPDX-License-Identifier: Unlicense */
-
-/* divide by three (based on routine from MPI and the GMP manual) */
-mp_err mp_div_3(const mp_int *a, mp_int *c, mp_digit *d)
-{
- mp_int q;
- mp_word w, t;
- mp_digit b;
- mp_err err;
- int ix;
-
- /* b = 2**MP_DIGIT_BIT / 3 */
- b = ((mp_word)1 << (mp_word)MP_DIGIT_BIT) / (mp_word)3;
-
- if ((err = mp_init_size(&q, a->used)) != MP_OKAY) {
- return err;
- }
-
- q.used = a->used;
- q.sign = a->sign;
- w = 0;
- for (ix = a->used - 1; ix >= 0; ix--) {
- w = (w << (mp_word)MP_DIGIT_BIT) | (mp_word)a->dp[ix];
-
- if (w >= 3u) {
- /* multiply w by [1/3] */
- t = (w * (mp_word)b) >> (mp_word)MP_DIGIT_BIT;
-
- /* now subtract 3 * [w/3] from w, to get the remainder */
- w -= t+t+t;
-
- /* fixup the remainder as required since
- * the optimization is not exact.
- */
- while (w >= 3u) {
- t += 1u;
- w -= 3u;
- }
- } else {
- t = 0;
- }
- q.dp[ix] = (mp_digit)t;
- }
-
- /* [optional] store the remainder */
- if (d != NULL) {
- *d = (mp_digit)w;
- }
-
- /* [optional] store the quotient */
- if (c != NULL) {
- mp_clamp(&q);
- mp_exch(&q, c);
- }
- mp_clear(&q);
-
- return err;
-}
-
-#endif
diff --git a/libtommath/bn_mp_expt_u32.c b/libtommath/bn_mp_expt_u32.c
deleted file mode 100644
index 2ab67ba..0000000
--- a/libtommath/bn_mp_expt_u32.c
+++ /dev/null
@@ -1,46 +0,0 @@
-#include "tommath_private.h"
-#ifdef BN_MP_EXPT_U32_C
-/* LibTomMath, multiple-precision integer library -- Tom St Denis */
-/* SPDX-License-Identifier: Unlicense */
-
-/* calculate c = a**b using a square-multiply algorithm */
-mp_err mp_expt_u32(const mp_int *a, uint32_t b, mp_int *c)
-{
- mp_err err;
-
- mp_int g;
-
- if ((err = mp_init_copy(&g, a)) != MP_OKAY) {
- return err;
- }
-
- /* set initial result */
- mp_set(c, 1uL);
-
- while (b > 0u) {
- /* if the bit is set multiply */
- if ((b & 1u) != 0u) {
- if ((err = mp_mul(c, &g, c)) != MP_OKAY) {
- goto LBL_ERR;
- }
- }
-
- /* square */
- if (b > 1u) {
- if ((err = mp_sqr(&g, &g)) != MP_OKAY) {
- goto LBL_ERR;
- }
- }
-
- /* shift to next bit */
- b >>= 1;
- }
-
- err = MP_OKAY;
-
-LBL_ERR:
- mp_clear(&g);
- return err;
-}
-
-#endif
diff --git a/libtommath/bn_mp_log_u32.c b/libtommath/bn_mp_log_u32.c
deleted file mode 100644
index b86d789..0000000
--- a/libtommath/bn_mp_log_u32.c
+++ /dev/null
@@ -1,180 +0,0 @@
-#include "tommath_private.h"
-#ifdef BN_MP_LOG_U32_C
-/* LibTomMath, multiple-precision integer library -- Tom St Denis */
-/* SPDX-License-Identifier: Unlicense */
-
-/* Compute log_{base}(a) */
-static mp_word s_pow(mp_word base, mp_word exponent)
-{
- mp_word result = 1u;
- while (exponent != 0u) {
- if ((exponent & 1u) == 1u) {
- result *= base;
- }
- exponent >>= 1;
- base *= base;
- }
-
- return result;
-}
-
-static mp_digit s_digit_ilogb(mp_digit base, mp_digit n)
-{
- mp_word bracket_low = 1u, bracket_mid, bracket_high, N;
- mp_digit ret, high = 1u, low = 0uL, mid;
-
- if (n < base) {
- return 0uL;
- }
- if (n == base) {
- return 1uL;
- }
-
- bracket_high = (mp_word) base ;
- N = (mp_word) n;
-
- while (bracket_high < N) {
- low = high;
- bracket_low = bracket_high;
- high <<= 1;
- bracket_high *= bracket_high;
- }
-
- while (((mp_digit)(high - low)) > 1u) {
- mid = (low + high) >> 1;
- bracket_mid = bracket_low * s_pow(base, (mp_word)(mid - low));
-
- if (N < bracket_mid) {
- high = mid ;
- bracket_high = bracket_mid ;
- }
- if (N > bracket_mid) {
- low = mid ;
- bracket_low = bracket_mid ;
- }
- if (N == bracket_mid) {
- return (mp_digit) mid;
- }
- }
-
- if (bracket_high == N) {
- ret = high;
- } else {
- ret = low;
- }
-
- return ret;
-}
-
-/* TODO: output could be "int" because the output of mp_radix_size is int, too,
- as is the output of mp_bitcount.
- With the same problem: max size is INT_MAX * MP_DIGIT not INT_MAX only!
-*/
-mp_err mp_log_u32(const mp_int *a, uint32_t base, uint32_t *c)
-{
- mp_err err;
- mp_ord cmp;
- uint32_t high, low, mid;
- mp_int bracket_low, bracket_high, bracket_mid, t, bi_base;
-
- err = MP_OKAY;
-
- if (a->sign == MP_NEG) {
- return MP_VAL;
- }
-
- if (MP_IS_ZERO(a)) {
- return MP_VAL;
- }
-
- if (base < 2u) {
- return MP_VAL;
- }
-
- /* A small shortcut for bases that are powers of two. */
- if ((base & (base - 1u)) == 0u) {
- int y, bit_count;
- for (y=0; (y < 7) && ((base & 1u) == 0u); y++) {
- base >>= 1;
- }
- bit_count = mp_count_bits(a) - 1;
- *c = (uint32_t)(bit_count/y);
- return MP_OKAY;
- }
-
- if (a->used == 1) {
- *c = (uint32_t)s_digit_ilogb(base, a->dp[0]);
- return err;
- }
-
- cmp = mp_cmp_d(a, base);
- if ((cmp == MP_LT) || (cmp == MP_EQ)) {
- *c = cmp == MP_EQ;
- return err;
- }
-
- if ((err =
- mp_init_multi(&bracket_low, &bracket_high,
- &bracket_mid, &t, &bi_base, NULL)) != MP_OKAY) {
- return err;
- }
-
- low = 0u;
- mp_set(&bracket_low, 1uL);
- high = 1u;
-
- mp_set(&bracket_high, base);
-
- /*
- A kind of Giant-step/baby-step algorithm.
- Idea shamelessly stolen from https://programmingpraxis.com/2010/05/07/integer-logarithms/2/
- The effect is asymptotic, hence needs benchmarks to test if the Giant-step should be skipped
- for small n.
- */
- while (mp_cmp(&bracket_high, a) == MP_LT) {
- low = high;
- if ((err = mp_copy(&bracket_high, &bracket_low)) != MP_OKAY) {
- goto LBL_ERR;
- }
- high <<= 1;
- if ((err = mp_sqr(&bracket_high, &bracket_high)) != MP_OKAY) {
- goto LBL_ERR;
- }
- }
- mp_set(&bi_base, base);
-
- while ((high - low) > 1u) {
- mid = (high + low) >> 1;
-
- if ((err = mp_expt_u32(&bi_base, (uint32_t)(mid - low), &t)) != MP_OKAY) {
- goto LBL_ERR;
- }
- if ((err = mp_mul(&bracket_low, &t, &bracket_mid)) != MP_OKAY) {
- goto LBL_ERR;
- }
- cmp = mp_cmp(a, &bracket_mid);
- if (cmp == MP_LT) {
- high = mid;
- mp_exch(&bracket_mid, &bracket_high);
- }
- if (cmp == MP_GT) {
- low = mid;
- mp_exch(&bracket_mid, &bracket_low);
- }
- if (cmp == MP_EQ) {
- *c = mid;
- goto LBL_END;
- }
- }
-
- *c = (mp_cmp(&bracket_high, a) == MP_EQ) ? high : low;
-
-LBL_END:
-LBL_ERR:
- mp_clear_multi(&bracket_low, &bracket_high, &bracket_mid,
- &t, &bi_base, NULL);
- return err;
-}
-
-
-#endif
diff --git a/libtommath/bn_mp_root_u32.c b/libtommath/bn_mp_root_u32.c
deleted file mode 100644
index ba65549..0000000
--- a/libtommath/bn_mp_root_u32.c
+++ /dev/null
@@ -1,139 +0,0 @@
-#include "tommath_private.h"
-#ifdef BN_MP_ROOT_U32_C
-/* LibTomMath, multiple-precision integer library -- Tom St Denis */
-/* SPDX-License-Identifier: Unlicense */
-
-/* find the n'th root of an integer
- *
- * Result found such that (c)**b <= a and (c+1)**b > a
- *
- * This algorithm uses Newton's approximation
- * x[i+1] = x[i] - f(x[i])/f'(x[i])
- * which will find the root in log(N) time where
- * each step involves a fair bit.
- */
-mp_err mp_root_u32(const mp_int *a, uint32_t b, mp_int *c)
-{
- mp_int t1, t2, t3, a_;
- mp_ord cmp;
- int ilog2;
- mp_err err;
-
- /* input must be positive if b is even */
- if (((b & 1u) == 0u) && (a->sign == MP_NEG)) {
- return MP_VAL;
- }
-
- if ((err = mp_init_multi(&t1, &t2, &t3, NULL)) != MP_OKAY) {
- return err;
- }
-
- /* if a is negative fudge the sign but keep track */
- a_ = *a;
- a_.sign = MP_ZPOS;
-
- /* Compute seed: 2^(log_2(n)/b + 2)*/
- ilog2 = mp_count_bits(a);
-
- /*
- If "b" is larger than INT_MAX it is also larger than
- log_2(n) because the bit-length of the "n" is measured
- with an int and hence the root is always < 2 (two).
- */
- if (b > (uint32_t)(INT_MAX/2)) {
- mp_set(c, 1uL);
- c->sign = a->sign;
- err = MP_OKAY;
- goto LBL_ERR;
- }
-
- /* "b" is smaller than INT_MAX, we can cast safely */
- if (ilog2 < (int)b) {
- mp_set(c, 1uL);
- c->sign = a->sign;
- err = MP_OKAY;
- goto LBL_ERR;
- }
- ilog2 = ilog2 / ((int)b);
- if (ilog2 == 0) {
- mp_set(c, 1uL);
- c->sign = a->sign;
- err = MP_OKAY;
- goto LBL_ERR;
- }
- /* Start value must be larger than root */
- ilog2 += 2;
- if ((err = mp_2expt(&t2,ilog2)) != MP_OKAY) goto LBL_ERR;
- do {
- /* t1 = t2 */
- if ((err = mp_copy(&t2, &t1)) != MP_OKAY) goto LBL_ERR;
-
- /* t2 = t1 - ((t1**b - a) / (b * t1**(b-1))) */
-
- /* t3 = t1**(b-1) */
- if ((err = mp_expt_u32(&t1, b - 1u, &t3)) != MP_OKAY) goto LBL_ERR;
-
- /* numerator */
- /* t2 = t1**b */
- if ((err = mp_mul(&t3, &t1, &t2)) != MP_OKAY) goto LBL_ERR;
-
- /* t2 = t1**b - a */
- if ((err = mp_sub(&t2, &a_, &t2)) != MP_OKAY) goto LBL_ERR;
-
- /* denominator */
- /* t3 = t1**(b-1) * b */
- if ((err = mp_mul_d(&t3, b, &t3)) != MP_OKAY) goto LBL_ERR;
-
- /* t3 = (t1**b - a)/(b * t1**(b-1)) */
- if ((err = mp_div(&t2, &t3, &t3, NULL)) != MP_OKAY) goto LBL_ERR;
-
- if ((err = mp_sub(&t1, &t3, &t2)) != MP_OKAY) goto LBL_ERR;
-
- /*
- Number of rounds is at most log_2(root). If it is more it
- got stuck, so break out of the loop and do the rest manually.
- */
- if (ilog2-- == 0) {
- break;
- }
- } while (mp_cmp(&t1, &t2) != MP_EQ);
-
- /* result can be off by a few so check */
- /* Loop beneath can overshoot by one if found root is smaller than actual root */
- for (;;) {
- if ((err = mp_expt_u32(&t1, b, &t2)) != MP_OKAY) goto LBL_ERR;
- cmp = mp_cmp(&t2, &a_);
- if (cmp == MP_EQ) {
- err = MP_OKAY;
- goto LBL_ERR;
- }
- if (cmp == MP_LT) {
- if ((err = mp_add_d(&t1, 1uL, &t1)) != MP_OKAY) goto LBL_ERR;
- } else {
- break;
- }
- }
- /* correct overshoot from above or from recurrence */
- for (;;) {
- if ((err = mp_expt_u32(&t1, b, &t2)) != MP_OKAY) goto LBL_ERR;
- if (mp_cmp(&t2, &a_) == MP_GT) {
- if ((err = mp_sub_d(&t1, 1uL, &t1)) != MP_OKAY) goto LBL_ERR;
- } else {
- break;
- }
- }
-
- /* set the result */
- mp_exch(&t1, c);
-
- /* set the sign of the result */
- c->sign = a->sign;
-
- err = MP_OKAY;
-
-LBL_ERR:
- mp_clear_multi(&t1, &t2, &t3, NULL);
- return err;
-}
-
-#endif
generated by cgit v0.12 at 2025-12-09 09:20:08 (GMT)