ref: 6ed61e66ac92b07c5424b0364d1c7589e1ddd437
parent: 689ec946d8650d617f72359e2c997ce2dfa61d3f
author: Jeff Bezanson <[email protected]>
date: Fri Feb 17 18:59:45 EST 2012
moving some flisp-specific code out of library
--- a/cvalues.c
+++ b/cvalues.c
@@ -1,3 +1,5 @@
+#include "operators.c"
+
#ifdef BITS64
#define NWORDS(sz) (((sz)+7)>>3)
#else
--- a/flisp.h
+++ b/flisp.h
@@ -208,6 +208,19 @@
} cvtable_t;
/* functions needed to implement the value interface (cvtable_t) */
+typedef enum { T_INT8, T_UINT8, T_INT16, T_UINT16, T_INT32, T_UINT32,+ T_INT64, T_UINT64, T_FLOAT, T_DOUBLE } numerictype_t;
+
+#define N_NUMTYPES ((int)T_DOUBLE+1)
+
+#ifdef BITS64
+# define T_LONG T_INT64
+# define T_ULONG T_UINT64
+#else
+# define T_LONG T_INT32
+# define T_ULONG T_UINT32
+#endif
+
value_t relocate_lispvalue(value_t v);
void print_traverse(value_t v);
void fl_print_chr(char c, ios_t *f);
@@ -350,6 +363,21 @@
value_t mk_wchar(int32_t n);
value_t return_from_uint64(uint64_t Uaccum);
value_t return_from_int64(int64_t Saccum);
+
+numerictype_t effective_numerictype(double r);
+double conv_to_double(void *data, numerictype_t tag);
+void conv_from_double(void *data, double d, numerictype_t tag);
+int64_t conv_to_int64(void *data, numerictype_t tag);
+uint64_t conv_to_uint64(void *data, numerictype_t tag);
+int32_t conv_to_int32(void *data, numerictype_t tag);
+uint32_t conv_to_uint32(void *data, numerictype_t tag);
+#ifdef BITS64
+#define conv_to_long conv_to_int64
+#define conv_to_ulong conv_to_uint64
+#else
+#define conv_to_long conv_to_int32
+#define conv_to_ulong conv_to_uint32
+#endif
typedef struct {char *name;
--- a/llt/Makefile
+++ b/llt/Makefile
@@ -1,14 +1,11 @@
CC = gcc
-SRCS = bitvector.c hashing.c socket.c timefuncs.c dblprint.c ptrhash.c \
- utf8.c ios.c operators.c cplxprint.c dirpath.c htable.c \
- bitvector-ops.c int2str.c dump.c random.c bswap.c memalign.c \
- swapreverse.c lltinit.c arraylist.c
+SRCS = bitvector.c hashing.c socket.c timefuncs.c ptrhash.c utf8.c ios.c \
+ dirpath.c htable.c bitvector-ops.c int2str.c dump.c random.c \
+ lltinit.c arraylist.c
OBJS = $(SRCS:%.c=%.o)
DOBJS = $(SRCS:%.c=%.do)
TARGET = libllt.a
-TESTSRC = unittest.c
-TESTER = llttest
# OS flags: LINUX, WIN32, MACOSX
# architecture flags: __CPU__=xxx, BITS64, ARCH_X86, ARCH_X86_64
@@ -34,12 +31,6 @@
rm -rf $(TARGET)
ar rs $(TARGET) $(OBJS)
-test:
- make clean
- make debug CFLAGS=-DENABLE_LLT_TEST
- gcc $(TESTSRC) $(TARGET) -o $(TESTER) -lm
- ./$(TESTER)
-
clean:
rm -f *.o
rm -f *.do
@@ -46,4 +37,3 @@
rm -f *~
rm -f core*
rm -f $(TARGET)
- rm -f $(TESTER)
--- a/llt/bswap.c
+++ /dev/null
@@ -1,87 +1,0 @@
-#include <stdlib.h>
-#include <string.h>
-#include <stddef.h>
-#include <alloca.h>
-#include "dtypes.h"
-#include "utils.h"
-
-void bswap_buffer(byte_t *data, size_t sz, size_t npts)
-{- size_t i, b;
- byte_t *el;
- byte_t temp;
-
- if (sz <= 1)
- return;
-
- switch (sz) {- case 8:
- for(i=0; i < npts; i++) {- ((u_int64_t*)data)[i] = bswap_64(((u_int64_t*)data)[i]);
- }
- break;
- case 4:
- for(i=0; i < npts; i++) {- ((u_int32_t*)data)[i] = bswap_32(((u_int32_t*)data)[i]);
- }
- break;
- case 2:
- for(i=0; i < npts; i++) {- ((u_int16_t*)data)[i] = bswap_16(((u_int16_t*)data)[i]);
- }
- break;
- default:
- for(i=0; i < sz * npts; i += sz) {- el = data + i;
- for(b=0; b < sz/2; b++) {- temp = el[b];
- el[b] = el[sz-b-1];
- el[sz-b-1] = temp;
- }
- }
- }
-}
-
-void bswap(byte_t *s, size_t n)
-{- unsigned int i;
- char temp;
-
- switch (n) {- case 8:
- *(u_int64_t*)s = bswap_64(*(u_int64_t*)s); break;
- case 4:
- *(u_int32_t*)s = bswap_32(*(u_int32_t*)s); break;
- case 2:
- *(u_int16_t*)s = bswap_16(*(u_int16_t*)s); break;
- case 1:
- break;
- default:
- for(i=0; i < n/2; i++) {- temp = s[i];
- s[i] = s[n-i-1];
- s[n-i-1] = temp;
- }
- }
-}
-
-void bswap_to(byte_t *dest, byte_t *src, size_t n)
-{- unsigned int i;
-
- switch (n) {- case 8:
- *(u_int64_t*)dest = bswap_64(*(u_int64_t*)src); break;
- case 4:
- *(u_int32_t*)dest = bswap_32(*(u_int32_t*)src); break;
- case 2:
- *(u_int16_t*)dest = bswap_16(*(u_int16_t*)src); break;
- case 1:
- break;
- default:
- for(i=0; i < n; i++) {- dest[i] = src[n-i-1];
- }
- }
-}
-
--- a/llt/cplxprint.c
+++ /dev/null
@@ -1,66 +1,0 @@
-#include <math.h>
-#include <stdlib.h>
-#include <string.h>
-#include "dtypes.h"
-#include "utils.h"
-
-void snprint_cplx(char *s, size_t cnt, double re, double im,
- // args to pass on to snprint_real
- int width, int dec,
- int max_digs_rt, int max_digs_lf,
- // print spaces around sign in a+bi
- int spflag)
-{- int fzr = (re==0) || rel_zero(re,im);
- int fzi = (im==0) || rel_zero(im,re);
- size_t len, sl;
- size_t space = cnt;
-
- s[0] = '\0';
- if (isnan(im) && fzr) {- if (space < 2) return;
- snprint_real(s, space-2, im, width, dec, max_digs_rt, max_digs_lf);
- strcat(s, "i");
- return;
- }
- if (!fzr || (fzr && fzi)) {- if (space < 4) return;
- snprint_real(s, space-4, re, width, dec, max_digs_rt, max_digs_lf);
- if ((im >= 0 || (isnan(im)&&!sign_bit(im))) && !fzi) {- if (spflag) {- strcat(s, " + ");
- }
- else {- strcat(s, "+");
- }
- }
- else if (!fzi) {- im = -im;
- if (spflag)
- strcat(s, " - ");
- else
- strcat(s, "-");
- }
- }
- if (!fzi) {- len = sl = strlen(s);
- if (im == -1) {- while ((long)(len-sl) < (long)(width-2) && len < (space-3))
- s[len++] = ' ';
- s[len] = '-';
- s[len+1] = 'i';
- s[len+2] = '\0';
- }
- else if (im == 1) {- while ((long)(len-sl) < (long)(width-1) && len < (space-2))
- s[len++] = ' ';
- s[len] = 'i';
- s[len+1] = '\0';
- }
- else {- snprint_real(s+len, space-len-2, im, width, dec,
- max_digs_rt, max_digs_lf);
- strcat(s, "i");
- }
- }
-}
--- a/llt/dblprint.c
+++ /dev/null
@@ -1,122 +1,0 @@
-#include <math.h>
-#include <stdlib.h>
-#include <string.h>
-#include <stdio.h>
-#include "dtypes.h"
-#include "ieee754.h"
-#include "utils.h"
-
-int double_exponent(double d)
-{- union ieee754_double dl;
-
- dl.d = d;
- return dl.ieee.exponent - IEEE754_DOUBLE_BIAS;
-}
-
-double double_mantissa(double d)
-{- union ieee754_double dl;
-
- dl.d = d;
- dl.ieee.exponent = IEEE754_DOUBLE_BIAS;
- dl.ieee.negative = 0;
- return dl.d;
-}
-
-int float_exponent(float f)
-{- union ieee754_float fl;
-
- fl.f = f;
- return fl.ieee.exponent - IEEE754_FLOAT_BIAS;
-}
-
-float float_mantissa(float f)
-{- union ieee754_float fl;
-
- fl.f = f;
- fl.ieee.exponent = IEEE754_FLOAT_BIAS;
- fl.ieee.negative = 0;
- return fl.f;
-}
-
-void snprint_real(char *s, size_t cnt, double r,
- int width, // printf field width, or 0
- int dec, // # decimal digits desired, recommend 16
- // # of zeros in .00...0x before using scientific notation
- // recommend 3-4 or so
- int max_digs_rt,
- // # of digits left of decimal before scientific notation
- // recommend 10
- int max_digs_lf)
-{- int mag;
- double fpart, temp;
- char format[8];
- char num_format[3];
- int sz, keepz=0;
-
- s[0] = '\0';
- if (width == -1) {- width = 0;
- keepz=1;
- }
- if (isnan(r)) {- if (sign_bit(r))
- strncpy(s, "-nan", cnt);
- else
- strncpy(s, "nan", cnt);
- return;
- }
- if (r == 0) {- strncpy(s, "0", cnt);
- return;
- }
-
- num_format[0] = 'l';
- num_format[2] = '\0';
-
- mag = double_exponent(r);
-
- mag = (int)(((double)mag)/LOG2_10 + 0.5);
- if (r == 0)
- mag = 0;
- if ((mag > max_digs_lf-1) || (mag < -max_digs_rt)) {- num_format[1] = 'e';
- temp = r/pow(10, mag); /* see if number will have a decimal */
- fpart = temp - floor(temp); /* when written in scientific notation */
- }
- else {- num_format[1] = 'f';
- fpart = r - floor(r);
- }
- if (fpart == 0)
- dec = 0;
- if (width == 0) {- snprintf(format, 8, "%%.%d%s", dec, num_format);
- }
- else {- snprintf(format, 8, "%%%d.%d%s", width, dec, num_format);
- }
- sz = snprintf(s, cnt, format, r);
- /* trim trailing zeros from fractions. not when using scientific
- notation, since we might have e.g. 1.2000e+100. also not when we
- need a specific output width */
- if (width == 0 && !keepz) {- if (sz > 2 && fpart && num_format[1]!='e') {- while (s[sz-1] == '0') {- s[sz-1]='\0';
- sz--;
- }
- // don't need trailing .
- if (s[sz-1] == '.') {- s[sz-1] = '\0';
- sz--;
- }
- }
- }
- // TODO. currently 1.1e20 prints as 1.1000000000000000e+20; be able to
- // get rid of all those zeros.
-}
--- a/llt/dtypes.h
+++ b/llt/dtypes.h
@@ -178,17 +178,4 @@
extern float F_PINF;
extern float F_NINF;
-typedef enum { T_INT8, T_UINT8, T_INT16, T_UINT16, T_INT32, T_UINT32,- T_INT64, T_UINT64, T_FLOAT, T_DOUBLE } numerictype_t;
-
-#define N_NUMTYPES ((int)T_DOUBLE+1)
-
-#ifdef BITS64
-# define T_LONG T_INT64
-# define T_ULONG T_UINT64
-#else
-# define T_LONG T_INT32
-# define T_ULONG T_UINT32
-#endif
-
#endif
--- a/llt/memalign.c
+++ /dev/null
@@ -1,50 +1,0 @@
-#include <stdlib.h>
-#include <string.h>
-#include <stddef.h>
-#include <alloca.h>
-#include "dtypes.h"
-#include "utils.h"
-
-#define ALIGNED_TO_ACTUAL(p) (((char*)p) - ((long*)p)[-1])
-
-static void *aligned_ptr(char *ptr, size_t align_size)
-{- char *ptr2, *aligned_ptr;
-
- ptr2 = ptr + sizeof(long);
- aligned_ptr = (char*)LLT_ALIGN(((uptrint_t)ptr2), align_size);
-
- ((long*)aligned_ptr)[-1] = (long)(aligned_ptr - ptr);
-
- return aligned_ptr;
-}
-
-/* align_size has to be a power of two */
-void *malloc_aligned(size_t size, size_t align_size)
-{- char *ptr;
-
- ptr = (char*)LLT_ALLOC(size + align_size-1 + sizeof(long));
- if (ptr == NULL)
- return NULL;
-
- return aligned_ptr(ptr, align_size);
-}
-
-void free_aligned(void *ptr)
-{- LLT_FREE(ALIGNED_TO_ACTUAL(ptr));
-}
-
-void *realloc_aligned(void *ptr, size_t size, size_t align_size)
-{- char *pnew;
-
- if (ptr != NULL)
- ptr = ALIGNED_TO_ACTUAL(ptr);
- pnew = LLT_REALLOC(ptr, size + align_size-1 + sizeof(long));
- if (pnew == NULL)
- return NULL;
-
- return aligned_ptr(pnew, align_size);
-}
--- a/llt/operators.c
+++ /dev/null
@@ -1,355 +1,0 @@
-#include <limits.h>
-#include <assert.h>
-#include "dtypes.h"
-#include "utils.h"
-#include "ieee754.h"
-
-// given a number, determine an appropriate type for storing it
-#if 0
-numerictype_t effective_numerictype(double r)
-{- double fp;
-
- fp = fpart(r);
- if (fp != 0 || r > U64_MAX || r < S64_MIN) {- if (r > FLT_MAX || r < -FLT_MAX || (fabs(r) < FLT_MIN)) {- return T_DOUBLE;
- }
- else {- return T_FLOAT;
- }
- }
- else if (r >= SCHAR_MIN && r <= SCHAR_MAX) {- return T_INT8;
- }
- else if (r >= SHRT_MIN && r <= SHRT_MAX) {- return T_INT16;
- }
- else if (r >= INT_MIN && r <= INT_MAX) {- return T_INT32;
- }
- else if (r <= S64_MAX) {- return T_INT64;
- }
- return T_UINT64;
-}
-#else
-// simpler version implementing a smaller preferred type repertoire
-numerictype_t effective_numerictype(double r)
-{- double fp;
-
- fp = fpart(r);
- if (fp != 0 || r > U64_MAX || r < S64_MIN) {- return T_DOUBLE;
- }
- else if (r >= INT_MIN && r <= INT_MAX) {- return T_INT32;
- }
- else if (r <= S64_MAX) {- return T_INT64;
- }
- return T_UINT64;
-}
-#endif
-
-double conv_to_double(void *data, numerictype_t tag)
-{- double d=0;
- switch (tag) {- case T_INT8: d = (double)*(int8_t*)data; break;
- case T_UINT8: d = (double)*(uint8_t*)data; break;
- case T_INT16: d = (double)*(int16_t*)data; break;
- case T_UINT16: d = (double)*(uint16_t*)data; break;
- case T_INT32: d = (double)*(int32_t*)data; break;
- case T_UINT32: d = (double)*(uint32_t*)data; break;
- case T_INT64:
- d = (double)*(int64_t*)data;
- if (d > 0 && *(int64_t*)data < 0) // can happen!
- d = -d;
- break;
- case T_UINT64: d = (double)*(uint64_t*)data; break;
- case T_FLOAT: d = (double)*(float*)data; break;
- case T_DOUBLE: return *(double*)data;
- }
- return d;
-}
-
-void conv_from_double(void *dest, double d, numerictype_t tag)
-{- switch (tag) {- case T_INT8: *(int8_t*)dest = d; break;
- case T_UINT8: *(uint8_t*)dest = d; break;
- case T_INT16: *(int16_t*)dest = d; break;
- case T_UINT16: *(uint16_t*)dest = d; break;
- case T_INT32: *(int32_t*)dest = d; break;
- case T_UINT32: *(uint32_t*)dest = d; break;
- case T_INT64:
- *(int64_t*)dest = d;
- if (d > 0 && *(int64_t*)dest < 0) // 0x8000000000000000 is a bitch
- *(int64_t*)dest = S64_MAX;
- break;
- case T_UINT64: *(uint64_t*)dest = (int64_t)d; break;
- case T_FLOAT: *(float*)dest = d; break;
- case T_DOUBLE: *(double*)dest = d; break;
- }
-}
-
-#define CONV_TO_INTTYPE(type) \
-type##_t conv_to_##type(void *data, numerictype_t tag) \
-{ \- type##_t i=0; \
- switch (tag) { \- case T_INT8: i = (type##_t)*(int8_t*)data; break; \
- case T_UINT8: i = (type##_t)*(uint8_t*)data; break; \
- case T_INT16: i = (type##_t)*(int16_t*)data; break; \
- case T_UINT16: i = (type##_t)*(uint16_t*)data; break; \
- case T_INT32: i = (type##_t)*(int32_t*)data; break; \
- case T_UINT32: i = (type##_t)*(uint32_t*)data; break; \
- case T_INT64: i = (type##_t)*(int64_t*)data; break; \
- case T_UINT64: i = (type##_t)*(uint64_t*)data; break; \
- case T_FLOAT: i = (type##_t)*(float*)data; break; \
- case T_DOUBLE: i = (type##_t)*(double*)data; break; \
- } \
- return i; \
-}
-
-CONV_TO_INTTYPE(int64)
-CONV_TO_INTTYPE(int32)
-CONV_TO_INTTYPE(uint32)
-
-// this is needed to work around a possible compiler bug
-// casting negative floats and doubles to uint64. you need
-// to cast to int64 first.
-uint64_t conv_to_uint64(void *data, numerictype_t tag)
-{- uint64_t i=0;
- switch (tag) {- case T_INT8: i = (uint64_t)*(int8_t*)data; break;
- case T_UINT8: i = (uint64_t)*(uint8_t*)data; break;
- case T_INT16: i = (uint64_t)*(int16_t*)data; break;
- case T_UINT16: i = (uint64_t)*(uint16_t*)data; break;
- case T_INT32: i = (uint64_t)*(int32_t*)data; break;
- case T_UINT32: i = (uint64_t)*(uint32_t*)data; break;
- case T_INT64: i = (uint64_t)*(int64_t*)data; break;
- case T_UINT64: i = (uint64_t)*(uint64_t*)data; break;
- case T_FLOAT:
- if (*(float*)data >= 0)
- i = (uint64_t)*(float*)data;
- else
- i = (uint64_t)(int64_t)*(float*)data;
- break;
- case T_DOUBLE:
- if (*(double*)data >= 0)
- i = (uint64_t)*(double*)data;
- else
- i = (uint64_t)(int64_t)*(double*)data;
- break;
- }
- return i;
-}
-
-int cmp_same_lt(void *a, void *b, numerictype_t tag)
-{- switch (tag) {- case T_INT8: return *(int8_t*)a < *(int8_t*)b;
- case T_UINT8: return *(uint8_t*)a < *(uint8_t*)b;
- case T_INT16: return *(int16_t*)a < *(int16_t*)b;
- case T_UINT16: return *(uint16_t*)a < *(uint16_t*)b;
- case T_INT32: return *(int32_t*)a < *(int32_t*)b;
- case T_UINT32: return *(uint32_t*)a < *(uint32_t*)b;
- case T_INT64: return *(int64_t*)a < *(int64_t*)b;
- case T_UINT64: return *(uint64_t*)a < *(uint64_t*)b;
- case T_FLOAT: return *(float*)a < *(float*)b;
- case T_DOUBLE: return *(double*)a < *(double*)b;
- }
- return 0;
-}
-
-int cmp_same_eq(void *a, void *b, numerictype_t tag)
-{- switch (tag) {- case T_INT8: return *(int8_t*)a == *(int8_t*)b;
- case T_UINT8: return *(uint8_t*)a == *(uint8_t*)b;
- case T_INT16: return *(int16_t*)a == *(int16_t*)b;
- case T_UINT16: return *(uint16_t*)a == *(uint16_t*)b;
- case T_INT32: return *(int32_t*)a == *(int32_t*)b;
- case T_UINT32: return *(uint32_t*)a == *(uint32_t*)b;
- case T_INT64: return *(int64_t*)a == *(int64_t*)b;
- case T_UINT64: return *(uint64_t*)a == *(uint64_t*)b;
- case T_FLOAT: return *(float*)a == *(float*)b;
- case T_DOUBLE: return *(double*)a == *(double*)b;
- }
- return 0;
-}
-
-int cmp_lt(void *a, numerictype_t atag, void *b, numerictype_t btag)
-{- if (atag==btag)
- return cmp_same_lt(a, b, atag);
-
- double da = conv_to_double(a, atag);
- double db = conv_to_double(b, btag);
-
- // casting to double will only get the wrong answer for big int64s
- // that differ in low bits
- if (da < db)
- return 1;
- if (db < da)
- return 0;
-
- if (atag == T_UINT64) {- // this is safe because if a had been bigger than S64_MAX,
- // we would already have concluded that it's bigger than b.
- if (btag == T_INT64) {- return ((int64_t)*(uint64_t*)a < *(int64_t*)b);
- }
- else if (btag == T_DOUBLE) {- if (db != db) return 0;
- return (*(uint64_t*)a < (uint64_t)*(double*)b);
- }
- }
- else if (atag == T_INT64) {- if (btag == T_UINT64) {- return (*(int64_t*)a < (int64_t)*(uint64_t*)b);
- }
- else if (btag == T_DOUBLE) {- if (db != db) return 0;
- return (*(int64_t*)a < (int64_t)*(double*)b);
- }
- }
- else if (btag == T_UINT64) {- if (atag == T_INT64) {- return ((int64_t)*(uint64_t*)b > *(int64_t*)a);
- }
- else if (atag == T_DOUBLE) {- if (da != da) return 0;
- return (*(uint64_t*)b > (uint64_t)*(double*)a);
- }
- }
- else if (btag == T_INT64) {- if (atag == T_UINT64) {- return (*(int64_t*)b > (int64_t)*(uint64_t*)a);
- }
- else if (atag == T_DOUBLE) {- if (da != da) return 0;
- return (*(int64_t*)b > (int64_t)*(double*)a);
- }
- }
- return 0;
-}
-
-int cmp_eq(void *a, numerictype_t atag, void *b, numerictype_t btag,
- int equalnans)
-{- if (atag==btag && (!equalnans || atag < T_FLOAT))
- return cmp_same_eq(a, b, atag);
-
- double da = conv_to_double(a, atag);
- double db = conv_to_double(b, btag);
-
- if ((int)atag >= T_FLOAT && (int)btag >= T_FLOAT) {- if (equalnans) {- return *(uint64_t*)&da == *(uint64_t*)&db;
- }
- return (da == db);
- }
-
- if (da != db)
- return 0;
-
- if (atag == T_UINT64) {- // this is safe because if a had been bigger than S64_MAX,
- // we would already have concluded that it's bigger than b.
- if (btag == T_INT64) {- return ((int64_t)*(uint64_t*)a == *(int64_t*)b);
- }
- else if (btag == T_DOUBLE) {- return (*(uint64_t*)a == (uint64_t)(int64_t)*(double*)b);
- }
- }
- else if (atag == T_INT64) {- if (btag == T_UINT64) {- return (*(int64_t*)a == (int64_t)*(uint64_t*)b);
- }
- else if (btag == T_DOUBLE) {- return (*(int64_t*)a == (int64_t)*(double*)b);
- }
- }
- else if (btag == T_UINT64) {- if (atag == T_INT64) {- return ((int64_t)*(uint64_t*)b == *(int64_t*)a);
- }
- else if (atag == T_DOUBLE) {- return (*(uint64_t*)b == (uint64_t)(int64_t)*(double*)a);
- }
- }
- else if (btag == T_INT64) {- if (atag == T_UINT64) {- return (*(int64_t*)b == (int64_t)*(uint64_t*)a);
- }
- else if (atag == T_DOUBLE) {- return (*(int64_t*)b == (int64_t)*(double*)a);
- }
- }
- return 1;
-}
-
-#ifdef ENABLE_LLT_TEST
-void test_operators()
-{- int8_t i8, i8b;
- uint8_t ui8, ui8b;
- int16_t i16, i16b;
- uint16_t ui16, ui16b;
- int32_t i32, i32b;
- uint32_t ui32, ui32b;
- int64_t i64, i64b;
- uint64_t ui64, ui64b;
- float f, fb;
- double d, db;
-
- ui64 = U64_MAX;
- ui64b = U64_MAX-1;
- i64 = S64_MIN;
- i64b = i64+1;
- d = (double)ui64;
- db = (double)i64b;
-
- assert(cmp_lt(&i64, T_INT64, &ui64, T_UINT64));
- assert(!cmp_lt(&ui64, T_UINT64, &i64, T_INT64));
- assert(cmp_lt(&i64, T_INT64, &ui64b, T_UINT64));
- assert(!cmp_lt(&ui64b, T_UINT64, &i64, T_INT64));
- assert(cmp_lt(&i64, T_INT64, &i64b, T_INT64));
- assert(!cmp_lt(&i64b, T_INT64, &i64, T_INT64));
-
- // try to compare a double too big to fit in an int64 with an
- // int64 requiring too much precision to fit in a double...
- // this case fails but it's very difficult/expensive to support
- //assert(cmp_lt(&ui64b, T_UINT64, &d, T_DOUBLE));
-
- i64 = S64_MAX;
- ui64 = S64_MAX-1;
- assert(cmp_lt(&ui64, T_UINT64, &i64, T_INT64));
- assert(!cmp_lt(&i64, T_INT64, &ui64, T_UINT64));
- i64 = S64_MAX-1;
- ui64 = S64_MAX;
- assert(cmp_lt(&i64, T_INT64, &ui64, T_UINT64));
- assert(!cmp_lt(&ui64, T_UINT64, &i64, T_INT64));
-
- d = DBL_MAXINT;
- i64 = DBL_MAXINT+100;
- assert(cmp_lt(&d, T_DOUBLE, &i64, T_INT64));
- assert(!cmp_lt(&i64, T_INT64, &d, T_DOUBLE));
- i64 = DBL_MAXINT+10;
- assert(cmp_lt(&d, T_DOUBLE, &i64, T_INT64));
- assert(!cmp_lt(&i64, T_INT64, &d, T_DOUBLE));
- i64 = DBL_MAXINT+1;
- assert(cmp_lt(&d, T_DOUBLE, &i64, T_INT64));
- assert(!cmp_lt(&i64, T_INT64, &d, T_DOUBLE));
-
- assert(!cmp_eq(&d, T_DOUBLE, &i64, T_INT64, 0));
- i64 = DBL_MAXINT;
- assert(cmp_eq(&d, T_DOUBLE, &i64, T_INT64, 0));
-}
-#endif
--- a/llt/swapreverse.c
+++ /dev/null
@@ -1,187 +1,0 @@
-#include <stdlib.h>
-#include <string.h>
-#include <stddef.h>
-#include <alloca.h>
-#include "dtypes.h"
-#include "utils.h"
-
-void memswap(char *a, char *b, size_t sz)
-{- int8_t i8;
- int32_t i32;
- int32_t *a4, *b4;
-
- if (sz < 4) {- while (sz--) {- i8 = *a;
- *a++ = *b;
- *b++ = i8;
- }
- }
- else {- while (sz & 0x3) {- i8 = *a;
- *a++ = *b;
- *b++ = i8;
- sz--;
- }
- a4 = (int32_t*)a;
- b4 = (int32_t*)b;
- sz >>= 2;
- while (sz--) {- i32 = *a4;
- *a4++ = *b4;
- *b4++ = i32;
- }
- }
-}
-
-void memreverse(char *a, size_t n, size_t elsz)
-{- int64_t i64, *pi64;
- int32_t i32, *pi32;
- int16_t i16, *pi16;
- int8_t i8;
- size_t i;
- char *temp;
- size_t eli, tot;
-
- if (n==0 || elsz==0) return;
- switch(elsz) {- case 16:
- pi64 = (int64_t*)a;
- for(i=0; i < n/2; i++) {- i64 = pi64[2*i];
- pi64[2*i] = pi64[2*(n-i-1)];
- pi64[2*(n-i-1)] = i64;
-
- i64 = pi64[2*i+1];
- pi64[2*i+1] = pi64[2*(n-i-1)+1];
- pi64[2*(n-i-1)+1] = i64;
- }
- break;
- case 8:
- pi64 = (int64_t*)a;
- for(i=0; i < n/2; i++) {- i64 = pi64[i];
- pi64[i] = pi64[n-i-1];
- pi64[n-i-1] = i64;
- }
- break;
- case 4:
- pi32 = (int32_t*)a;
- for(i=0; i < n/2; i++) {- i32 = pi32[i];
- pi32[i] = pi32[n-i-1];
- pi32[n-i-1] = i32;
- }
- break;
- case 2:
- pi16 = (int16_t*)a;
- for(i=0; i < n/2; i++) {- i16 = pi16[i];
- pi16[i] = pi16[n-i-1];
- pi16[n-i-1] = i16;
- }
- break;
- case 1:
- for(i=0; i < n/2; i++) {- i8 = a[i];
- a[i] = a[n-i-1];
- a[n-i-1] = i8;
- }
- break;
- default:
- tot = n*elsz;
- if (elsz < 4097)
- temp = alloca(elsz);
- else
- temp = malloc(elsz);
-
- if (temp != NULL) {- for(i=0, eli=0; i < n/2; i++, eli+=elsz) {- memcpy(temp, &a[eli], elsz);
- memcpy(&a[eli], &a[tot-eli-elsz], elsz);
- memcpy(&a[tot-eli-elsz], temp, elsz);
- }
-
- if (elsz >= 4097)
- free(temp);
- }
- break;
- }
-}
-
-void memreverse_to(char *dest, char *a, size_t n, size_t elsz)
-{- int64_t *pi64, *di64;
- int32_t *pi32, *di32;
- int16_t *pi16, *di16;
- size_t i;
- size_t eli, tot;
- if (n==0 || elsz==0) return;
- switch(elsz) {- case 16:
- pi64 = (int64_t*)a;
- di64 = (int64_t*)dest;
- for(i=0; i < n/2; i++) {- di64[2*i] = pi64[2*(n-i-1)];
- di64[2*(n-i-1)] = pi64[2*i];
-
- di64[2*i+1] = pi64[2*(n-i-1)+1];
- di64[2*(n-i-1)+1] = pi64[2*i+1];
- }
- if (n&0x1) {- di64[2*i] = pi64[2*i];
- di64[2*i+1] = pi64[2*i+1];
- }
- break;
- case 8:
- pi64 = (int64_t*)a;
- di64 = (int64_t*)dest;
- for(i=0; i < n/2; i++) {- di64[i] = pi64[n-i-1];
- di64[n-i-1] = pi64[i];
- }
- if (n&0x1)
- di64[i] = pi64[i];
- break;
- case 4:
- pi32 = (int32_t*)a;
- di32 = (int32_t*)dest;
- for(i=0; i < n/2; i++) {- di32[i] = pi32[n-i-1];
- di32[n-i-1] = pi32[i];
- }
- if (n&0x1)
- di32[i] = pi32[i];
- break;
- case 2:
- pi16 = (int16_t*)a;
- di16 = (int16_t*)dest;
- for(i=0; i < n/2; i++) {- di16[i] = pi16[n-i-1];
- di16[n-i-1] = pi16[i];
- }
- if (n&0x1)
- di16[i] = pi16[i];
- break;
- case 1:
- for(i=0; i < n/2; i++) {- dest[i] = a[n-i-1];
- dest[n-i-1] = a[i];
- }
- if (n&0x1)
- dest[i] = a[i];
- break;
- default:
- tot = n*elsz;
- for(i=0, eli=0; i < n/2; i++, eli+=elsz) {- memcpy(&dest[eli], &a[tot - eli - elsz], elsz);
- memcpy(&dest[tot - eli - elsz], &a[eli], elsz);
- }
- if (n&0x1)
- memcpy(&dest[eli], &a[eli], elsz);
- break;
- }
-}
--- a/llt/unittest.c
+++ /dev/null
@@ -1,142 +1,0 @@
-#include <stdlib.h>
-#include <stdarg.h>
-#include <stdio.h>
-#include <wchar.h>
-#include "llt.h"
-
-int main()
-{- llt_init();
-
- test_dblprint();
- test_operators();
-
- /*
- char *buf = malloc(20000);
- char *buf2 = malloc(20000);
- FILE *f = fopen("textread.m","rb");- int i=0;
- while (!feof(f))
- buf[i++] = fgetc(f);
- buf[i-1] = '\0';
- int len = i-1;
- double t0 = clock_now();
- int j=0;
- for(i=0; i < 20000; i++) {- //j+=u8_charnum(buf,len);
- u8_reverse(buf2, buf, len);
- }
- printf("textread took %.4f sec (%d)\n", clock_now()-t0, j);-
- FILE *f2 = fopen("u8.txt","rb");- i=0;
- while (!feof(f2))
- buf[i++] = fgetc(f2);
- buf[i-1] = '\0';
- len = i-1;
- t0 = clock_now();
- j=0;
- for(i=0; i < 20000; i++) {- //j+=u8_charnum(buf,len);
- u8_reverse(buf2, buf, len);
- }
- printf("u8 took %.4f sec (%d)\n\n", clock_now()-t0, j);- */
-
- test_ios();
-
- return 0;
-}
-
-static void prettycplx(double r, double i)
-{- char str[64];
- snprint_cplx(str, sizeof(str), r, i, 0, 16, 3, 10, 1);
- fputs(str, stdout);
- fputc('\n', stdout);-}
-
-static void prettyreal(double r)
-{- char str[64];
- snprint_real(str, sizeof(str), r, 0, 16, 3, 10);
- fputs(str, stdout);
- fputc('\n', stdout);-}
-
-void test_dblprint()
-{- char str[64];
-
- prettycplx(0,0);
- prettycplx(1,0);
- prettycplx(0,1);
- prettycplx(1,1);
- prettycplx(-1,0);
- prettycplx(0,-1);
- prettycplx(1,-1);
- prettycplx(-1,1);
- prettycplx(-1,-1);
- prettycplx(2,0);
- prettycplx(0,2);
- prettycplx(2,2);
- prettycplx(-2,0);
- prettycplx(0,-2);
- prettycplx(2,-2);
- prettycplx(-2,2);
- prettycplx(-2,-2);
-
- prettyreal(1.5);
- prettyreal(1.1);
- prettyreal(1.1e-100);
- prettyreal(1.1e20);
- prettyreal(123456789);
- prettyreal(1234567890);
- prettyreal(12345678901);
- prettyreal(-12345678901);
- prettyreal(12345678901223);
- prettyreal(-12345678901223);
- prettyreal(.02);
- prettyreal(.002);
- prettyreal(.0002);
- prettyreal(-.0002);
- prettyreal(.00002);
- prettyreal(-.00002);
-
- prettyreal(1.0/0);
- prettyreal(-1.0/0);
- prettyreal(strtod("nan",NULL));- prettyreal(0.0/0);
- prettyreal(-0.0/0);
-
- prettyreal(DBL_EPSILON);
-}
-
-void test_ios()
-{- ios_t *out = ios_stdout();
- ios_t *in = ios_stdin();
-
- ios_putc('a', out);- ios_putc('b', out);- ios_putc('\n', out);-
- char c[80];
- size_t i=0;
- ios_t sts;
- ios_str(&sts, "Test string.");
- c[i++] = ios_getc(&sts);
- c[i++] = ios_getc(&sts);
- c[i++] = ios_getc(&sts);
- c[i++] = '\0';
- printf("got: \"%s\"\n", c);-
- ios_t ms;
- ios_mem(&ms, 10);
- int j;
- for(j=0; j < 16; j++)
- ios_puts("passersby were amazed by the ", &ms);- size_t bs;
- char *bigstr = ios_takebuf(&ms, &bs);
- printf("got: \"%s\" (size %d)\n", bigstr, bs);-}
--- a/llt/unittest.c.1
+++ /dev/null
@@ -1,81 +1,0 @@
-#include <stdlib.h>
-#include <stdarg.h>
-#include <stdio.h>
-#include <wchar.h>
-#include "llt.h"
-
-int main()
-{- llt_init();
-
- test_dblprint();
- test_operators();
-
- return 0;
-}
-
-static void prettycplx(double r, double i)
-{- char str[64];
- snprint_cplx(str, sizeof(str), r, i, 0, 16, 3, 10, 1);
- fputs(str, stdout);
- fputc('\n', stdout);-}
-
-static void prettyreal(double r)
-{- char str[64];
- snprint_real(str, sizeof(str), r, 0, 16, 3, 10);
- fputs(str, stdout);
- fputc('\n', stdout);-}
-
-void test_dblprint()
-{- char str[64];
-
- dbl_tolerance(1e-12);
-
- prettycplx(0,0);
- prettycplx(1,0);
- prettycplx(0,1);
- prettycplx(1,1);
- prettycplx(-1,0);
- prettycplx(0,-1);
- prettycplx(1,-1);
- prettycplx(-1,1);
- prettycplx(-1,-1);
- prettycplx(2,0);
- prettycplx(0,2);
- prettycplx(2,2);
- prettycplx(-2,0);
- prettycplx(0,-2);
- prettycplx(2,-2);
- prettycplx(-2,2);
- prettycplx(-2,-2);
-
- prettyreal(1.5);
- prettyreal(1.1);
- prettyreal(1.1e-100);
- prettyreal(1.1e20);
- prettyreal(123456789);
- prettyreal(1234567890);
- prettyreal(12345678901);
- prettyreal(-12345678901);
- prettyreal(12345678901223);
- prettyreal(-12345678901223);
- prettyreal(.02);
- prettyreal(.002);
- prettyreal(.0002);
- prettyreal(-.0002);
- prettyreal(.00002);
- prettyreal(-.00002);
-
- prettyreal(1.0/0);
- prettyreal(-1.0/0);
- prettyreal(strtod("nan",NULL));- prettyreal(0.0/0);
- prettyreal(-0.0/0);
-
- prettyreal(DBL_EPSILON);
-}
--- a/llt/utils.h
+++ b/llt/utils.h
@@ -1,78 +1,9 @@
#ifndef __UTILS_H_
#define __UTILS_H_
-/* these functions byteswap any-size units --------------------- */
-void bswap(byte_t *s, size_t n);
-void bswap_to(byte_t *dest, byte_t *src, size_t n);
-void bswap_buffer(byte_t *data, size_t sz, size_t npts);
-/* ------------------------------------------------------------- */
-
-/* reverse the order of elements of any size, in place or out of place */
-/* n is the number of elements */
-void memreverse(char *a, size_t n, size_t elsz);
-void memreverse_to(char *dest, char *a, size_t n, size_t elsz);
-
-/* swap the contents of two buffers */
-void memswap(char *a, char *b, size_t sz);
-
-/* allocating aligned blocks ----------------------------------- */
-void *malloc_aligned(size_t size, size_t align_size);
-void free_aligned(void *ptr);
-void *realloc_aligned(void *ptr, size_t size, size_t align_size);
-/* ------------------------------------------------------------- */
-
-int double_exponent(double d);
-double double_mantissa(double d);
-int float_exponent(float f);
-float float_mantissa(float f);
-void snprint_real(char *s, size_t cnt, double r,
- int width, // printf field width, or 0
- int dec, // # decimal digits desired, recommend 16
- // # of zeros in .00...0x before using scientific notation
- // recommend 3-4 or so
- int max_digs_rt,
- // # of digits left of decimal before scientific notation
- // recommend 10
- int max_digs_lf);
-void snprint_cplx(char *s, size_t cnt, double re, double im,
- // args to pass on to snprint_real
- int width, int dec,
- int max_digs_rt, int max_digs_lf,
- // print spaces around sign in a+bi
- int spflag);
-
char *uint2str(char *dest, size_t len, uint64_t num, uint32_t base);
int str2int(char *str, size_t len, int64_t *res, uint32_t base);
int isdigit_base(char c, int base);
-
-extern double trunc(double x);
-
-STATIC_INLINE double fpart(double arg)
-{- return arg - trunc(arg);
-}
-
-#define ipart(x) trunc(x)
-
-numerictype_t effective_numerictype(double r);
-double conv_to_double(void *data, numerictype_t tag);
-void conv_from_double(void *data, double d, numerictype_t tag);
-int64_t conv_to_int64(void *data, numerictype_t tag);
-uint64_t conv_to_uint64(void *data, numerictype_t tag);
-int32_t conv_to_int32(void *data, numerictype_t tag);
-uint32_t conv_to_uint32(void *data, numerictype_t tag);
-#ifdef BITS64
-#define conv_to_long conv_to_int64
-#define conv_to_ulong conv_to_uint64
-#else
-#define conv_to_long conv_to_int32
-#define conv_to_ulong conv_to_uint32
-#endif
-int cmp_same_lt(void *a, void *b, numerictype_t tag);
-int cmp_same_eq(void *a, void *b, numerictype_t tag);
-int cmp_lt(void *a, numerictype_t atag, void *b, numerictype_t btag);
-int cmp_eq(void *a, numerictype_t atag, void *b, numerictype_t btag,
- int equalnans);
#ifdef ARCH_X86_64
# define LEGACY_REGS "=Q"
--- /dev/null
+++ b/operators.c
@@ -1,0 +1,362 @@
+#include <limits.h>
+#include <assert.h>
+#include "dtypes.h"
+#include "utils.h"
+#include "ieee754.h"
+
+extern double trunc(double x);
+
+STATIC_INLINE double fpart(double arg)
+{+ return arg - trunc(arg);
+}
+
+// given a number, determine an appropriate type for storing it
+#if 0
+numerictype_t effective_numerictype(double r)
+{+ double fp;
+
+ fp = fpart(r);
+ if (fp != 0 || r > U64_MAX || r < S64_MIN) {+ if (r > FLT_MAX || r < -FLT_MAX || (fabs(r) < FLT_MIN)) {+ return T_DOUBLE;
+ }
+ else {+ return T_FLOAT;
+ }
+ }
+ else if (r >= SCHAR_MIN && r <= SCHAR_MAX) {+ return T_INT8;
+ }
+ else if (r >= SHRT_MIN && r <= SHRT_MAX) {+ return T_INT16;
+ }
+ else if (r >= INT_MIN && r <= INT_MAX) {+ return T_INT32;
+ }
+ else if (r <= S64_MAX) {+ return T_INT64;
+ }
+ return T_UINT64;
+}
+#else
+// simpler version implementing a smaller preferred type repertoire
+numerictype_t effective_numerictype(double r)
+{+ double fp;
+
+ fp = fpart(r);
+ if (fp != 0 || r > U64_MAX || r < S64_MIN) {+ return T_DOUBLE;
+ }
+ else if (r >= INT_MIN && r <= INT_MAX) {+ return T_INT32;
+ }
+ else if (r <= S64_MAX) {+ return T_INT64;
+ }
+ return T_UINT64;
+}
+#endif
+
+double conv_to_double(void *data, numerictype_t tag)
+{+ double d=0;
+ switch (tag) {+ case T_INT8: d = (double)*(int8_t*)data; break;
+ case T_UINT8: d = (double)*(uint8_t*)data; break;
+ case T_INT16: d = (double)*(int16_t*)data; break;
+ case T_UINT16: d = (double)*(uint16_t*)data; break;
+ case T_INT32: d = (double)*(int32_t*)data; break;
+ case T_UINT32: d = (double)*(uint32_t*)data; break;
+ case T_INT64:
+ d = (double)*(int64_t*)data;
+ if (d > 0 && *(int64_t*)data < 0) // can happen!
+ d = -d;
+ break;
+ case T_UINT64: d = (double)*(uint64_t*)data; break;
+ case T_FLOAT: d = (double)*(float*)data; break;
+ case T_DOUBLE: return *(double*)data;
+ }
+ return d;
+}
+
+void conv_from_double(void *dest, double d, numerictype_t tag)
+{+ switch (tag) {+ case T_INT8: *(int8_t*)dest = d; break;
+ case T_UINT8: *(uint8_t*)dest = d; break;
+ case T_INT16: *(int16_t*)dest = d; break;
+ case T_UINT16: *(uint16_t*)dest = d; break;
+ case T_INT32: *(int32_t*)dest = d; break;
+ case T_UINT32: *(uint32_t*)dest = d; break;
+ case T_INT64:
+ *(int64_t*)dest = d;
+ if (d > 0 && *(int64_t*)dest < 0) // 0x8000000000000000 is a bitch
+ *(int64_t*)dest = S64_MAX;
+ break;
+ case T_UINT64: *(uint64_t*)dest = (int64_t)d; break;
+ case T_FLOAT: *(float*)dest = d; break;
+ case T_DOUBLE: *(double*)dest = d; break;
+ }
+}
+
+#define CONV_TO_INTTYPE(type) \
+type##_t conv_to_##type(void *data, numerictype_t tag) \
+{ \+ type##_t i=0; \
+ switch (tag) { \+ case T_INT8: i = (type##_t)*(int8_t*)data; break; \
+ case T_UINT8: i = (type##_t)*(uint8_t*)data; break; \
+ case T_INT16: i = (type##_t)*(int16_t*)data; break; \
+ case T_UINT16: i = (type##_t)*(uint16_t*)data; break; \
+ case T_INT32: i = (type##_t)*(int32_t*)data; break; \
+ case T_UINT32: i = (type##_t)*(uint32_t*)data; break; \
+ case T_INT64: i = (type##_t)*(int64_t*)data; break; \
+ case T_UINT64: i = (type##_t)*(uint64_t*)data; break; \
+ case T_FLOAT: i = (type##_t)*(float*)data; break; \
+ case T_DOUBLE: i = (type##_t)*(double*)data; break; \
+ } \
+ return i; \
+}
+
+CONV_TO_INTTYPE(int64)
+CONV_TO_INTTYPE(int32)
+CONV_TO_INTTYPE(uint32)
+
+// this is needed to work around a possible compiler bug
+// casting negative floats and doubles to uint64. you need
+// to cast to int64 first.
+uint64_t conv_to_uint64(void *data, numerictype_t tag)
+{+ uint64_t i=0;
+ switch (tag) {+ case T_INT8: i = (uint64_t)*(int8_t*)data; break;
+ case T_UINT8: i = (uint64_t)*(uint8_t*)data; break;
+ case T_INT16: i = (uint64_t)*(int16_t*)data; break;
+ case T_UINT16: i = (uint64_t)*(uint16_t*)data; break;
+ case T_INT32: i = (uint64_t)*(int32_t*)data; break;
+ case T_UINT32: i = (uint64_t)*(uint32_t*)data; break;
+ case T_INT64: i = (uint64_t)*(int64_t*)data; break;
+ case T_UINT64: i = (uint64_t)*(uint64_t*)data; break;
+ case T_FLOAT:
+ if (*(float*)data >= 0)
+ i = (uint64_t)*(float*)data;
+ else
+ i = (uint64_t)(int64_t)*(float*)data;
+ break;
+ case T_DOUBLE:
+ if (*(double*)data >= 0)
+ i = (uint64_t)*(double*)data;
+ else
+ i = (uint64_t)(int64_t)*(double*)data;
+ break;
+ }
+ return i;
+}
+
+int cmp_same_lt(void *a, void *b, numerictype_t tag)
+{+ switch (tag) {+ case T_INT8: return *(int8_t*)a < *(int8_t*)b;
+ case T_UINT8: return *(uint8_t*)a < *(uint8_t*)b;
+ case T_INT16: return *(int16_t*)a < *(int16_t*)b;
+ case T_UINT16: return *(uint16_t*)a < *(uint16_t*)b;
+ case T_INT32: return *(int32_t*)a < *(int32_t*)b;
+ case T_UINT32: return *(uint32_t*)a < *(uint32_t*)b;
+ case T_INT64: return *(int64_t*)a < *(int64_t*)b;
+ case T_UINT64: return *(uint64_t*)a < *(uint64_t*)b;
+ case T_FLOAT: return *(float*)a < *(float*)b;
+ case T_DOUBLE: return *(double*)a < *(double*)b;
+ }
+ return 0;
+}
+
+int cmp_same_eq(void *a, void *b, numerictype_t tag)
+{+ switch (tag) {+ case T_INT8: return *(int8_t*)a == *(int8_t*)b;
+ case T_UINT8: return *(uint8_t*)a == *(uint8_t*)b;
+ case T_INT16: return *(int16_t*)a == *(int16_t*)b;
+ case T_UINT16: return *(uint16_t*)a == *(uint16_t*)b;
+ case T_INT32: return *(int32_t*)a == *(int32_t*)b;
+ case T_UINT32: return *(uint32_t*)a == *(uint32_t*)b;
+ case T_INT64: return *(int64_t*)a == *(int64_t*)b;
+ case T_UINT64: return *(uint64_t*)a == *(uint64_t*)b;
+ case T_FLOAT: return *(float*)a == *(float*)b;
+ case T_DOUBLE: return *(double*)a == *(double*)b;
+ }
+ return 0;
+}
+
+int cmp_lt(void *a, numerictype_t atag, void *b, numerictype_t btag)
+{+ if (atag==btag)
+ return cmp_same_lt(a, b, atag);
+
+ double da = conv_to_double(a, atag);
+ double db = conv_to_double(b, btag);
+
+ // casting to double will only get the wrong answer for big int64s
+ // that differ in low bits
+ if (da < db)
+ return 1;
+ if (db < da)
+ return 0;
+
+ if (atag == T_UINT64) {+ // this is safe because if a had been bigger than S64_MAX,
+ // we would already have concluded that it's bigger than b.
+ if (btag == T_INT64) {+ return ((int64_t)*(uint64_t*)a < *(int64_t*)b);
+ }
+ else if (btag == T_DOUBLE) {+ if (db != db) return 0;
+ return (*(uint64_t*)a < (uint64_t)*(double*)b);
+ }
+ }
+ else if (atag == T_INT64) {+ if (btag == T_UINT64) {+ return (*(int64_t*)a < (int64_t)*(uint64_t*)b);
+ }
+ else if (btag == T_DOUBLE) {+ if (db != db) return 0;
+ return (*(int64_t*)a < (int64_t)*(double*)b);
+ }
+ }
+ else if (btag == T_UINT64) {+ if (atag == T_INT64) {+ return ((int64_t)*(uint64_t*)b > *(int64_t*)a);
+ }
+ else if (atag == T_DOUBLE) {+ if (da != da) return 0;
+ return (*(uint64_t*)b > (uint64_t)*(double*)a);
+ }
+ }
+ else if (btag == T_INT64) {+ if (atag == T_UINT64) {+ return (*(int64_t*)b > (int64_t)*(uint64_t*)a);
+ }
+ else if (atag == T_DOUBLE) {+ if (da != da) return 0;
+ return (*(int64_t*)b > (int64_t)*(double*)a);
+ }
+ }
+ return 0;
+}
+
+int cmp_eq(void *a, numerictype_t atag, void *b, numerictype_t btag,
+ int equalnans)
+{+ if (atag==btag && (!equalnans || atag < T_FLOAT))
+ return cmp_same_eq(a, b, atag);
+
+ double da = conv_to_double(a, atag);
+ double db = conv_to_double(b, btag);
+
+ if ((int)atag >= T_FLOAT && (int)btag >= T_FLOAT) {+ if (equalnans) {+ return *(uint64_t*)&da == *(uint64_t*)&db;
+ }
+ return (da == db);
+ }
+
+ if (da != db)
+ return 0;
+
+ if (atag == T_UINT64) {+ // this is safe because if a had been bigger than S64_MAX,
+ // we would already have concluded that it's bigger than b.
+ if (btag == T_INT64) {+ return ((int64_t)*(uint64_t*)a == *(int64_t*)b);
+ }
+ else if (btag == T_DOUBLE) {+ return (*(uint64_t*)a == (uint64_t)(int64_t)*(double*)b);
+ }
+ }
+ else if (atag == T_INT64) {+ if (btag == T_UINT64) {+ return (*(int64_t*)a == (int64_t)*(uint64_t*)b);
+ }
+ else if (btag == T_DOUBLE) {+ return (*(int64_t*)a == (int64_t)*(double*)b);
+ }
+ }
+ else if (btag == T_UINT64) {+ if (atag == T_INT64) {+ return ((int64_t)*(uint64_t*)b == *(int64_t*)a);
+ }
+ else if (atag == T_DOUBLE) {+ return (*(uint64_t*)b == (uint64_t)(int64_t)*(double*)a);
+ }
+ }
+ else if (btag == T_INT64) {+ if (atag == T_UINT64) {+ return (*(int64_t*)b == (int64_t)*(uint64_t*)a);
+ }
+ else if (atag == T_DOUBLE) {+ return (*(int64_t*)b == (int64_t)*(double*)a);
+ }
+ }
+ return 1;
+}
+
+#ifdef ENABLE_LLT_TEST
+void test_operators()
+{+ int8_t i8, i8b;
+ uint8_t ui8, ui8b;
+ int16_t i16, i16b;
+ uint16_t ui16, ui16b;
+ int32_t i32, i32b;
+ uint32_t ui32, ui32b;
+ int64_t i64, i64b;
+ uint64_t ui64, ui64b;
+ float f, fb;
+ double d, db;
+
+ ui64 = U64_MAX;
+ ui64b = U64_MAX-1;
+ i64 = S64_MIN;
+ i64b = i64+1;
+ d = (double)ui64;
+ db = (double)i64b;
+
+ assert(cmp_lt(&i64, T_INT64, &ui64, T_UINT64));
+ assert(!cmp_lt(&ui64, T_UINT64, &i64, T_INT64));
+ assert(cmp_lt(&i64, T_INT64, &ui64b, T_UINT64));
+ assert(!cmp_lt(&ui64b, T_UINT64, &i64, T_INT64));
+ assert(cmp_lt(&i64, T_INT64, &i64b, T_INT64));
+ assert(!cmp_lt(&i64b, T_INT64, &i64, T_INT64));
+
+ // try to compare a double too big to fit in an int64 with an
+ // int64 requiring too much precision to fit in a double...
+ // this case fails but it's very difficult/expensive to support
+ //assert(cmp_lt(&ui64b, T_UINT64, &d, T_DOUBLE));
+
+ i64 = S64_MAX;
+ ui64 = S64_MAX-1;
+ assert(cmp_lt(&ui64, T_UINT64, &i64, T_INT64));
+ assert(!cmp_lt(&i64, T_INT64, &ui64, T_UINT64));
+ i64 = S64_MAX-1;
+ ui64 = S64_MAX;
+ assert(cmp_lt(&i64, T_INT64, &ui64, T_UINT64));
+ assert(!cmp_lt(&ui64, T_UINT64, &i64, T_INT64));
+
+ d = DBL_MAXINT;
+ i64 = DBL_MAXINT+100;
+ assert(cmp_lt(&d, T_DOUBLE, &i64, T_INT64));
+ assert(!cmp_lt(&i64, T_INT64, &d, T_DOUBLE));
+ i64 = DBL_MAXINT+10;
+ assert(cmp_lt(&d, T_DOUBLE, &i64, T_INT64));
+ assert(!cmp_lt(&i64, T_INT64, &d, T_DOUBLE));
+ i64 = DBL_MAXINT+1;
+ assert(cmp_lt(&d, T_DOUBLE, &i64, T_INT64));
+ assert(!cmp_lt(&i64, T_INT64, &d, T_DOUBLE));
+
+ assert(!cmp_eq(&d, T_DOUBLE, &i64, T_INT64, 0));
+ i64 = DBL_MAXINT;
+ assert(cmp_eq(&d, T_DOUBLE, &i64, T_INT64, 0));
+}
+#endif
--- a/print.c
+++ b/print.c
@@ -1,3 +1,5 @@
+#include "ieee754.h"
+
extern void *memrchr(const void *s, int c, size_t n);
static htable_t printconses;
@@ -516,6 +518,93 @@
}
}
outc('"', f);+}
+
+int double_exponent(double d)
+{+ union ieee754_double dl;
+
+ dl.d = d;
+ return dl.ieee.exponent - IEEE754_DOUBLE_BIAS;
+}
+
+void snprint_real(char *s, size_t cnt, double r,
+ int width, // printf field width, or 0
+ int dec, // # decimal digits desired, recommend 16
+ // # of zeros in .00...0x before using scientific notation
+ // recommend 3-4 or so
+ int max_digs_rt,
+ // # of digits left of decimal before scientific notation
+ // recommend 10
+ int max_digs_lf)
+{+ int mag;
+ double fpart, temp;
+ char format[8];
+ char num_format[3];
+ int sz, keepz=0;
+
+ s[0] = '\0';
+ if (width == -1) {+ width = 0;
+ keepz=1;
+ }
+ if (isnan(r)) {+ if (sign_bit(r))
+ strncpy(s, "-nan", cnt);
+ else
+ strncpy(s, "nan", cnt);
+ return;
+ }
+ if (r == 0) {+ strncpy(s, "0", cnt);
+ return;
+ }
+
+ num_format[0] = 'l';
+ num_format[2] = '\0';
+
+ mag = double_exponent(r);
+
+ mag = (int)(((double)mag)/LOG2_10 + 0.5);
+ if (r == 0)
+ mag = 0;
+ if ((mag > max_digs_lf-1) || (mag < -max_digs_rt)) {+ num_format[1] = 'e';
+ temp = r/pow(10, mag); /* see if number will have a decimal */
+ fpart = temp - floor(temp); /* when written in scientific notation */
+ }
+ else {+ num_format[1] = 'f';
+ fpart = r - floor(r);
+ }
+ if (fpart == 0)
+ dec = 0;
+ if (width == 0) {+ snprintf(format, 8, "%%.%d%s", dec, num_format);
+ }
+ else {+ snprintf(format, 8, "%%%d.%d%s", width, dec, num_format);
+ }
+ sz = snprintf(s, cnt, format, r);
+ /* trim trailing zeros from fractions. not when using scientific
+ notation, since we might have e.g. 1.2000e+100. also not when we
+ need a specific output width */
+ if (width == 0 && !keepz) {+ if (sz > 2 && fpart && num_format[1]!='e') {+ while (s[sz-1] == '0') {+ s[sz-1]='\0';
+ sz--;
+ }
+ // don't need trailing .
+ if (s[sz-1] == '.') {+ s[sz-1] = '\0';
+ sz--;
+ }
+ }
+ }
+ // TODO. currently 1.1e20 prints as 1.1000000000000000e+20; be able to
+ // get rid of all those zeros.
}
static numerictype_t sym_to_numtype(value_t type);