libzahl

big integer library
git clone git://git.suckless.org/libzahl
Log | Files | Refs | README | LICENSE

commit 69787408e645ca3c28c13bf2de3e86f632f5549e
parent a41bcd346b1a5624c2c600d995a23bf701da118f
Author: Mattias Andrée <maandree@kth.se>
Date:   Fri, 22 Apr 2016 06:46:50 +0200

Add support for benchmark against TomsFastMath

Signed-off-by: Mattias Andrée <maandree@kth.se>

Diffstat:
MMakefile | 8+++++++-
MREADME | 7+++++--
Abench/libtfm.h | 293+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
3 files changed, 305 insertions(+), 3 deletions(-)

diff --git a/Makefile b/Makefile @@ -81,6 +81,8 @@ BENCHMARK_LIB_tommath = -ltommath BENCHMARK_LIB_libtommath = -ltommath BENCHMARK_LIB_gmp = -lgmp BENCHMARK_LIB_libgmp = -lgmp +BENCHMARK_LIB_tfm = libtfm.a +BENCHMARK_LIB_libtfm = libtfm.a BENCHMARK_DEP_ = libzahl.a BENCHMARK_DEP_zahl = libzahl.a @@ -89,11 +91,15 @@ BENCHMARK_DEP_tommath = bench/libtommath.h BENCHMARK_DEP_libtommath = bench/libtommath.h BENCHMARK_DEP_gmp = bench/libgmp.h BENCHMARK_DEP_libgmp = bench/libgmp.h +BENCHMARK_DEP_tfm = bench/libtfm.h +BENCHMARK_DEP_libtfm = bench/libtfm.h BENCHMARK_CPP_tommath = '-DBENCHMARK_LIB="libtommath.h"' BENCHMARK_CPP_libtommath = '-DBENCHMARK_LIB="libtommath.h"' BENCHMARK_CPP_gmp = '-DBENCHMARK_LIB="libgmp.h"' BENCHMARK_CPP_libgmp = '-DBENCHMARK_LIB="libgmp.h"' +BENCHMARK_CPP_tfm = '-DBENCHMARK_LIB="libtfm.h"' +BENCHMARK_CPP_libtfm = '-DBENCHMARK_LIB="libtfm.h"' CPPFLAGS += $(BENCHMARK_CPP_$(BENCHMARK_LIB)) @@ -145,6 +151,6 @@ uninstall: -cd -- "$(DESTDIR)$(MANPREFIX)/man7" && rm $(MAN7) clean: - -rm -- *.o *.su *.a *.so test test-random.c 2>/dev/null + -rm -- *.o *.su *.a *.so test test-random.c config.mk 2>/dev/null .PHONY: all check clean install uninstall diff --git a/README b/README @@ -26,8 +26,11 @@ DESCRIPTION RATIONALE GMP MP cannot be used for rubust programs. LibTomMath is too slow, probably because of all memory allocations, - and has an nonintuitive API. Hebimath is promising, but - I think it can be done better. + and has a nonintuitive API. TomsFastMath has an a + nonintuitive API, has limited precision (selected at + compile-time), and has limited functionality. All the + above are also bloated. Hebimath is promising, but I + think it can be done better. NOTES libzahl is currently not thread-safe. diff --git a/bench/libtfm.h b/bench/libtfm.h @@ -0,0 +1,293 @@ +#include <tfm.h> + +#include <stddef.h> +#include <setjmp.h> +#include <stdint.h> +#include <stdio.h> + +typedef fp_int z_t[1]; + +static z_t _0, _1, _a, _b; +static int _tmp; + +static void +fp_set_int(z_t a, uint32_t d) +{ + a->dp[0] = d; + a->used = !!d; + a->sign = 0; +} + +static void +zsetup(jmp_buf env) +{ + (void) env; + fp_init(_0); + fp_init(_1); + fp_init(_a); + fp_init(_b); + fp_set_int(_0, 0); + fp_set_int(_1, 1); +} + +static void +zunsetup(void) +{ + /* Do nothing */ +} + +#define FAST_RANDOM 0 +#define SECURE_RANDOM 0 +#define DEFAULT_RANDOM 0 +#define FASTEST_RANDOM 0 +#define LIBC_RAND_RANDOM 0 +#define LIBC_RANDOM_RANDOM 0 +#define LIBC_RAND48_RANDOM 0 +#define QUASIUNIFORM 0 +#define UNIFORM 1 +#define MODUNIFORM 2 + +#define zperror(x) ((void)0) +#define zinit(a) fp_init(a) +#define zfree(a) ((void)a) + +#define zset(r, a) fp_copy(a, r) +#define zadd_unsigned(r, a, b) (zabs(_a, a), zabs(_b, b), fp_add(_a, _b, r)) +#define zsub_unsigned(r, a, b) (zabs(_a, a), zabs(_b, b), fp_sub(_a, _b, r)) +#define zadd(r, a, b) fp_add(a, b, r) +#define zsub(r, a, b) fp_sub(a, b, r) +#define zeven(a) fp_iseven(a) +#define zodd(a) fp_isodd(a) +#define zeven_nonzero(a) fp_iseven(a) +#define zodd_nonzero(a) fp_isodd(a) +#define zzero(a) fp_iszero(a) +#define zsignum(a) fp_cmp(a, _0) +#define zbits(a) fp_count_bits(a) +#define zlsb(a) fp_cnt_lsb(a) +#define zlsh(r, a, b) fp_mul_2d(a, b, r) +#define zrsh(r, a, b) fp_div_2d(a, b, r, 0) +#define ztrunc(r, a, b) fp_mod_2d(a, b, r) +#define zcmpmag(a, b) fp_cmp_mag(a, b) +#define zcmp(a, b) fp_cmp(a, b) +#define zcmpi(a, b) (zseti(_b, b), fp_cmp(a, _b)) +#define zcmpu(a, b) (zsetu(_b, b), fp_cmp(a, _b)) +#define zgcd(r, a, b) fp_gcd(a, b, r) +#define zmul(r, a, b) fp_mul(a, b, r) +#define zsqr(r, a) fp_sqr(a, r) +#define zmodmul(r, a, b, m) fp_mulmod(a, b, m, r) +#define zmodsqr(r, a, m) fp_sqrmod(a, m, r) +#define zpow(r, a, b) zpowu(r, a, b->used ? b->dp[0] : 0) +#define zmodpow(r, a, b, m) fp_exptmod(a, b, m, r) +#define zmodpowu(r, a, b, m) (fp_set_int(_b, b), fp_exptmod(a, _b, m, r)) +#define zsets(a, s) fp_read_radix(a, s, 10) +#define zstr_length(a, b) (fp_radix_size(a, b, &_tmp), _tmp) +#define zstr(a, s) fp_toradix(a, s, 10) +#define zptest(w, a, t) fp_isprime_ex(a, t) /* Note, the witness is not returned. */ +#define zload(a, s) fp_read_signed_bin(a, (unsigned char *)s, _tmp) +#define zdiv(r, a, b) fp_div(a, b, r, 0) +#define zmod(r, a, b) fp_mod(a, b, r) +#define zdivmod(q, r, a, b) fp_div(a, b, q, r) + +static inline void +zneg(z_t r, z_t a) +{ + fp_neg(a, r); +} + +static inline void +zabs(z_t r, z_t a) +{ + fp_abs(a, r); +} + +static void +zswap(z_t a, z_t b) +{ + z_t t; + fp_copy(a, t); + fp_copy(b, a); + fp_copy(t, b); +} + +static void +zand(z_t r, z_t a, z_t b) +{ + int i; + for (i = 0; i < a->used && i < b->used; i++) + r->dp[i] = a->dp[i] & b->dp[i]; + r->used = i; + while (r->used && !r->dp[r->used]) + r->used--; + r->sign = a->sign & b->sign; +} + +static void +zor(z_t r, z_t a, z_t b) +{ + int i; + for (i = 0; i < a->used && i < b->used; i++) + r->dp[i] = a->dp[i] | b->dp[i]; + for (; i < a->used; i++) + r->dp[i] = a->dp[i]; + for (; i < b->used; i++) + r->dp[i] = b->dp[i]; + r->used = i; + r->sign = a->sign | b->sign; +} + +static void +zxor(z_t r, z_t a, z_t b) +{ + int i; + for (i = 0; i < a->used && i < b->used; i++) + r->dp[i] = a->dp[i] ^ b->dp[i]; + for (; i < a->used; i++) + r->dp[i] = a->dp[i]; + for (; i < b->used; i++) + r->dp[i] = b->dp[i]; + r->used = i; + while (r->used && !r->dp[r->used]) + r->used--; + r->sign = a->sign ^ b->sign; +} + +static int +zsave(z_t a, char *buf) +{ + _tmp = buf ? fp_signed_bin_size(a) : (fp_to_signed_bin(a, (unsigned char *)buf), 0); + return _tmp; +} + +static void +zsetu(z_t r, unsigned long long int val) +{ + uint32_t high = (uint32_t)(val >> 32); + uint32_t low = (uint32_t)val; + + if (high) { + fp_set_int(r, high); + fp_set_int(_a, low); + fp_lshd(r, 32); + zadd(r, r, _a); + } else { + fp_set_int(r, low); + } + +} + +static void +zseti(z_t r, long long int val) +{ + if (val >= 0) { + zsetu(r, (unsigned long long int)val); + } else { + zsetu(r, (unsigned long long int)-val); + zneg(r, r); + } +} + +static void +zsplit(z_t high, z_t low, z_t a, size_t brk) +{ + if (low == a) { + zrsh(high, a, brk); + ztrunc(low, a, brk); + } else { + ztrunc(low, a, brk); + zrsh(high, a, brk); + } +} + +static void +znot(z_t r, z_t a) +{ + fp_2expt(_a, (int)zbits(a)); + fp_sub_d(_a, 1, _a); + zand(r, a, _a); + zneg(r, r); +} + +static int +zbtest(z_t a, size_t bit) +{ + fp_2expt(_b, (int)bit); + zand(_b, a, _b); + return !zzero(_b); +} + +static void +zbset(z_t r, z_t a, size_t bit, int mode) +{ + if (mode > 0) { + fp_2expt(_b, (int)bit); + zor(r, a, _b); + } else if (mode < 0 || zbtest(a, bit)) { + fp_2expt(_b, (int)bit); + zxor(r, a, _b); + } +} + +static void +zrand(z_t r, int dev, int dist, z_t n) +{ + static int gave_up = 0; + int bits; + (void) dev; + + if (zzero(n)) { + fp_zero(r); + return; + } + if (zsignum(n) < 0) { + return; + } + + bits = zbits(n); + + switch (dist) { + case QUASIUNIFORM: + fp_rand(r, bits); + zadd(r, r, _1); + zmul(r, r, n); + zrsh(r, r, bits); + break; + + case UNIFORM: + if (!gave_up) { + gave_up = 1; + printf("I'm sorry, this is too difficult, I give up.\n"); + } + break; + + case MODUNIFORM: + fp_rand(r, bits); + if (zcmp(r, n) > 0) + zsub(r, r, n); + break; + + default: + abort(); + } +} + +static void +zpowu(z_t r, z_t a, unsigned long long int b) +{ + z_t product, factor; + int neg = zsignum(a) < 0; + + zinit(product); + zinit(factor); + zsetu(product, 1); + zset(factor, a); + + for (; b; b >>= 1) { + if (b & 1) + zmul(product, product, factor); + zsqr(factor, factor); + } + + zset(r, product); + if (neg) + zneg(r, r); +}