#endif
#ifdef HAVE_INTTYPES_H
-#include <inttypes.h>
+# include <inttypes.h>
#endif
-
-#include <string.h>
-#include <stdlib.h>
-#include <stdio.h>
-#include <assert.h>
-#ifdef HAVE_ALLOCA_H
-# include <alloca.h>
+#ifdef HAVE_STRING_H
+# include <string.h>
#endif
-#ifdef HAVE_MALLOC_H
-# include <malloc.h>
+#ifdef HAVE_STDLIB_H
+# include <stdlib.h>
#endif
+#include <stdio.h>
+#include <assert.h>
+
+#include "xmalloc.h"
typedef uint32_t UINT32;
#endif
#endif
-
-/********
- * globals
- ********/
+/**
+ * possible float states
+ */
typedef enum {
- NORMAL,
- ZERO,
- SUBNORMAL,
- INF,
- NAN,
+ NORMAL, /**< normal representation, implicit 1 */
+ ZERO, /**< +/-0 */
+ SUBNORMAL, /**< denormals, implicit 0 */
+ INF, /**< +/-oo */
+ NAN, /**< Not A Number */
} value_class_t;
+/** A descriptor for an IEEE float value. */
typedef struct {
- char exponent_size;
- char mantissa_size;
- value_class_t class;
+ unsigned char exponent_size; /**< size of exponent in bits */
+ unsigned char mantissa_size; /**< size of mantissa in bits */
+ value_class_t clss; /**< state of this float */
} descriptor_t;
-#define CLEAR_BUFFER(buffer) memset(buffer, 0, CALC_BUFFER_SIZE)
+#define CLEAR_BUFFER(buffer) memset(buffer, 0, calc_buffer_size)
/* because variable sized structs are impossible, the internal
* value is represented as a pseudo-struct char array, addressed
* by macros
* struct {
* char sign; // 0 for positive, 1 for negative
- * char exp[VALUE_SIZE];
- * char mant[VALUE_SIZE];
+ * char exp[value_size];
+ * char mant[value_size];
* descriptor_t desc;
* };
*/
#define _sign(a) (((char*)a)[SIGN_POS])
#define _exp(a) (&((char*)a)[EXPONENT_POS])
#define _mant(a) (&((char*)a)[MANTISSA_POS])
-#define _desc(a) (*(descriptor_t *)&((char*)a)[DESCRIPTOR_POS])
+#define _desc(a) (*(descriptor_t *)&((char *)a)[DESCRIPTOR_POS])
-#define _save_result(x) memcpy((x), sc_get_buffer(), VALUE_SIZE)
-#define _shift_right(x, y, b) sc_shr((x), (y), VALUE_SIZE*4, 0, (b))
-#define _shift_left(x, y, b) sc_shl((x), (y), VALUE_SIZE*4, 0, (b))
+#define _save_result(x) memcpy((x), sc_get_buffer(), value_size)
+#define _shift_right(x, y, b) sc_shr((x), (y), value_size*4, 0, (b))
+#define _shift_left(x, y, b) sc_shl((x), (y), value_size*4, 0, (b))
-#define FC_DEFINE1(code) char* fc_##code(const void *a, void *result) \
- { \
- return _calc((const char*)a, NULL, FC_##code, (char*)result); \
- }
+#define FC_DEFINE1(code) \
+char *fc_##code(const void *a, void *result) { \
+ return _calc((const char*)a, NULL, FC_##code, (char*)result); \
+}
-#define FC_DEFINE2(code) char* fc_##code(const void *a, const void *b, void *result) \
- { \
- return _calc((const char*)a, (const char*)b, FC_##code, (char*)result); \
- }
+#define FC_DEFINE2(code) \
+char *fc_##code(const void *a, const void *b, void *result) { \
+ return _calc((const char*)a, (const char*)b, FC_##code, (char*)result); \
+}
#define FUNC_PTR(code) fc_##code
static char *calc_buffer = NULL;
-static fc_rounding_mode_t ROUNDING_MODE;
+static fc_rounding_mode_t rounding_mode;
-static int CALC_BUFFER_SIZE;
-static int VALUE_SIZE;
+static int calc_buffer_size;
+static int value_size;
static int SIGN_POS;
static int EXPONENT_POS;
static int MANTISSA_POS;
}
#endif
-/* pack machine-like */
+/** pack machine-like */
static char* _pack(const char *int_float, char *packed)
{
char *shift_val;
char *temp;
char *val_buffer;
- temp = alloca(VALUE_SIZE);
- shift_val = alloca(VALUE_SIZE);
+ temp = alloca(value_size);
+ shift_val = alloca(value_size);
- switch (_desc(int_float).class) {
+ switch (_desc(int_float).clss) {
case NAN:
- val_buffer = alloca(CALC_BUFFER_SIZE);
+ val_buffer = alloca(calc_buffer_size);
fc_get_qnan(_desc(int_float).exponent_size, _desc(int_float).mantissa_size, val_buffer);
int_float = val_buffer;
break;
case INF:
- val_buffer = alloca(CALC_BUFFER_SIZE);
+ val_buffer = alloca(calc_buffer_size);
fc_get_plusinf(_desc(int_float).exponent_size, _desc(int_float).mantissa_size, val_buffer);
_sign(val_buffer) = _sign(int_float);
int_float = val_buffer;
char lsb, guard, round, round_dir = 0;
char *temp;
- temp = alloca(VALUE_SIZE);
+ temp = alloca(value_size);
/* +2: save two rounding bits at the end */
hsb = 2 + _desc(in_val).mantissa_size - sc_get_highest_set_bit(_mant(in_val)) - 1;
memcpy(&_desc(out_val), &_desc(in_val), sizeof(descriptor_t));
}
- _desc(out_val).class = NORMAL;
+ _desc(out_val).clss = NORMAL;
- /* mantissa all zeroes, so zero exponent (because of explicit one)*/
+ /* mantissa all zeros, so zero exponent (because of explicit one)*/
if (hsb == 2 + _desc(in_val).mantissa_size)
{
sc_val_from_ulong(0, _exp(out_val));
/* denormalized means exponent of zero */
sc_val_from_ulong(0, _exp(out_val));
- _desc(out_val).class = SUBNORMAL;
+ _desc(out_val).clss = SUBNORMAL;
}
/* perform rounding by adding a value that clears the guard bit and the round bit
guard = (lsb&0x2)>>1;
round = lsb&0x1;
- switch (ROUNDING_MODE)
+ switch (rounding_mode)
{
case FC_TONEAREST:
/* round to nearest representable value, if in doubt choose the version
}
/* could have rounded down to zero */
- if (sc_is_zero(_mant(out_val)) && (_desc(out_val).class == SUBNORMAL))
- _desc(out_val).class = ZERO;
+ if (sc_is_zero(_mant(out_val)) && (_desc(out_val).clss == SUBNORMAL))
+ _desc(out_val).clss = ZERO;
/* check for rounding overflow */
hsb = 2 + _desc(out_val).mantissa_size - sc_get_highest_set_bit(_mant(out_val)) - 1;
- if ((_desc(out_val).class != SUBNORMAL) && (hsb < -1))
+ if ((_desc(out_val).clss != SUBNORMAL) && (hsb < -1))
{
sc_val_from_ulong(1, temp);
_shift_right(_mant(out_val), temp, _mant(out_val));
sc_add(_exp(out_val), temp, _exp(out_val));
}
- else if ((_desc(out_val).class == SUBNORMAL) && (hsb == -1))
+ else if ((_desc(out_val).clss == SUBNORMAL) && (hsb == -1))
{
- /* overflow caused the matissa to be normal again,
+ /* overflow caused the mantissa to be normal again,
* so adapt the exponent accordingly */
sc_val_from_ulong(1, temp);
sc_add(_exp(out_val), temp, _exp(out_val));
- _desc(out_val).class = NORMAL;
+ _desc(out_val).clss = NORMAL;
}
/* no further rounding is needed, because rounding overflow means
* the carry of the original rounding was propagated all the way
if (_sign(out_val) == 0)
{
/* value is positive */
- switch (ROUNDING_MODE) {
+ switch (rounding_mode) {
case FC_TONEAREST:
case FC_TOPOSITIVE:
- _desc(out_val).class = INF;
+ _desc(out_val).clss = INF;
break;
case FC_TONEGATIVE:
}
} else {
/* value is negative */
- switch (ROUNDING_MODE) {
+ switch (rounding_mode) {
case FC_TONEAREST:
case FC_TONEGATIVE:
- _desc(out_val).class = INF;
+ _desc(out_val).clss = INF;
break;
case FC_TOPOSITIVE:
return out_val;
}
-/*
+/**
+ * Operations involving NaN's must return NaN
+ */
+#define handle_NAN(a, b, result) \
+do { \
+ if (_desc(a).clss == NAN) { \
+ if (a != result) memcpy(result, a, calc_buffer_size); \
+ return result; \
+ } \
+ if (_desc(b).clss == NAN) { \
+ if (b != result) memcpy(result, b, calc_buffer_size); \
+ return result; \
+ } \
+}while (0)
+
+
+/**
* calculate a + b, where a is the value with the bigger exponent
*/
-static char* _add(const char* a, const char* b, char* result)
+static char* _fadd(const char* a, const char* b, char* result)
{
char *temp;
char *exp_diff;
char sign;
char sticky;
- if (_desc(a).class == NAN) {
- if (a != result) memcpy(result, a, CALC_BUFFER_SIZE);
- return result;
- }
- if (_desc(b).class == NAN) {
- if (b != result) memcpy(result, b, CALC_BUFFER_SIZE);
- return result;
- }
+ handle_NAN(a, b, result);
/* make sure result has a descriptor */
if (result != a && result != b)
/* determine if this is an addition or subtraction */
sign = _sign(a) ^ _sign(b);
- /* produce nan on inf - inf */
- if (sign && (_desc(a).class == INF) && (_desc(b).class == INF))
+ /* produce NaN on inf - inf */
+ if (sign && (_desc(a).clss == INF) && (_desc(b).clss == INF))
return fc_get_qnan(_desc(a).exponent_size, _desc(b).mantissa_size, result);
- temp = alloca(VALUE_SIZE);
- exp_diff = alloca(VALUE_SIZE);
+ temp = alloca(value_size);
+ exp_diff = alloca(value_size);
/* get exponent difference */
sc_sub(_exp(a), _exp(b), exp_diff);
if (sign && sc_val_to_long(exp_diff) == 0) {
switch (sc_comp(_mant(a), _mant(b))) {
case 1: /* a > b */
- if (_sign(a)) _sign(result) = 1; /* abs(a) is bigger and a is negative */
- else _sign(result) = 0;
+ _sign(result) = _sign(a); /* abs(a) is bigger and a is negative */
break;
case 0: /* a == b */
- if (ROUNDING_MODE == FC_TONEGATIVE)
- _sign(result) = 1;
- else
- _sign(result) = 0;
+ _sign(result) = (rounding_mode == FC_TONEGATIVE);
break;
case -1: /* a < b */
- if (_sign(b)) _sign(result) = 1; /* abs(b) is bigger and b is negative */
- else _sign(result) = 0;
+ _sign(result) = _sign(b); /* abs(b) is bigger and b is negative */
break;
default:
/* can't be reached */
break;
}
- } else {
- _sign(result) = _sign(a);
}
+ else
+ _sign(result) = _sign(a);
/* sign has been taken care of, check for special cases */
- if (_desc(a).class == ZERO) {
- if (b != result) memcpy(result+SIGN_POS+1, b+SIGN_POS+1, CALC_BUFFER_SIZE-SIGN_POS-1);
+ if (_desc(a).clss == ZERO) {
+ if (b != result) memcpy(result+SIGN_POS+1, b+SIGN_POS+1, calc_buffer_size-SIGN_POS-1);
return result;
}
- if (_desc(b).class == ZERO) {
- if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, CALC_BUFFER_SIZE-SIGN_POS-1);
+ if (_desc(b).clss == ZERO) {
+ if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, calc_buffer_size-SIGN_POS-1);
return result;
}
- if (_desc(a).class == INF) {
- if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, CALC_BUFFER_SIZE-SIGN_POS-1);
+ if (_desc(a).clss == INF) {
+ if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, calc_buffer_size-SIGN_POS-1);
return result;
}
- if (_desc(b).class == INF) {
- if (b != result) memcpy(result+SIGN_POS+1, b+SIGN_POS+1, CALC_BUFFER_SIZE-SIGN_POS-1);
+ if (_desc(b).clss == INF) {
+ if (b != result) memcpy(result+SIGN_POS+1, b+SIGN_POS+1, calc_buffer_size-SIGN_POS-1);
return result;
}
/* shift the smaller value to the right to align the radix point */
/* subnormals have their radix point shifted to the right,
* take care of this first */
- if ((_desc(b).class == SUBNORMAL) && (_desc(a).class != SUBNORMAL))
+ if ((_desc(b).clss == SUBNORMAL) && (_desc(a).clss != SUBNORMAL))
{
sc_val_from_ulong(1, temp);
sc_sub(exp_diff, temp, exp_diff);
/* if subtracting a little more than the represented value or adding a little
* more than the represented value to a negative value this, in addition to the
* still set sticky bit, takes account of the 'little more' */
- char *temp1 = alloca(CALC_BUFFER_SIZE);
+ char *temp1 = alloca(calc_buffer_size);
sc_val_from_ulong(1, temp1);
sc_add(temp, temp1, temp);
}
/* _normalize expects a 'normal' radix point, adding two subnormals
* results in a subnormal radix point -> shifting before normalizing */
- if ((_desc(a).class == SUBNORMAL) && (_desc(b).class == SUBNORMAL))
+ if ((_desc(a).clss == SUBNORMAL) && (_desc(b).clss == SUBNORMAL))
{
sc_val_from_ulong(1, NULL);
_shift_left(_mant(result), sc_get_buffer(), _mant(result));
}
/* resulting exponent is the bigger one */
- memmove(_exp(result), _exp(a), VALUE_SIZE);
+ memmove(_exp(result), _exp(a), value_size);
return _normalize(result, result, sticky);
}
-static char* _mul(const char* a, const char* b, char* result)
+/**
+ * calculate a * b
+ */
+static char* _fmul(const char* a, const char* b, char* result)
{
char *temp;
- if (_desc(a).class == NAN) {
- if (a != result) memcpy(result, a, CALC_BUFFER_SIZE);
- return result;
- }
- if (_desc(b).class == NAN) {
- if (b != result) memcpy(result, b, CALC_BUFFER_SIZE);
- return result;
- }
+ handle_NAN(a, b, result);
- temp = alloca(VALUE_SIZE);
+ temp = alloca(value_size);
if (result != a && result != b)
memcpy(&_desc(result), &_desc(a), sizeof(descriptor_t));
_sign(result) = _sign(a) ^ _sign(b);
- /* produce nan on 0 * inf */
- if (_desc(a).class == ZERO) {
- if (_desc(b).class == INF)
+ /* produce NaN on 0 * inf */
+ if (_desc(a).clss == ZERO) {
+ if (_desc(b).clss == INF)
fc_get_qnan(_desc(a).exponent_size, _desc(a).mantissa_size, result);
else
- if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, CALC_BUFFER_SIZE-1);
+ if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, calc_buffer_size-1);
return result;
}
- if (_desc(b).class == ZERO) {
- if (_desc(a).class == INF)
+ if (_desc(b).clss == ZERO) {
+ if (_desc(a).clss == INF)
fc_get_qnan(_desc(a).exponent_size, _desc(a).mantissa_size, result);
else
- if (b != result) memcpy(result+SIGN_POS+1, b+SIGN_POS+1, CALC_BUFFER_SIZE-1);
+ if (b != result) memcpy(result+SIGN_POS+1, b+SIGN_POS+1, calc_buffer_size-1);
return result;
}
- if (_desc(a).class == INF) {
- if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, CALC_BUFFER_SIZE-1);
+ if (_desc(a).clss == INF) {
+ if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, calc_buffer_size-1);
return result;
}
- if (_desc(b).class == INF) {
- if (b != result) memcpy(result+SIGN_POS+1, b+SIGN_POS+1, CALC_BUFFER_SIZE-1);
+ if (_desc(b).clss == INF) {
+ if (b != result) memcpy(result+SIGN_POS+1, b+SIGN_POS+1, calc_buffer_size-1);
return result;
}
sc_sub(_exp(result), temp, _exp(result));
/* mixed normal, subnormal values introduce an error of 1, correct it */
- if ((_desc(a).class == SUBNORMAL) ^ (_desc(b).class == SUBNORMAL))
+ if ((_desc(a).clss == SUBNORMAL) ^ (_desc(b).clss == SUBNORMAL))
{
sc_val_from_ulong(1, temp);
sc_add(_exp(result), temp, _exp(result));
return _normalize(result, result, sc_had_carry());
}
-static char* _div(const char* a, const char* b, char* result)
+/**
+ * calculate a / b
+ */
+static char* _fdiv(const char* a, const char* b, char* result)
{
char *temp, *dividend;
- if (_desc(a).class == NAN) {
- if (a != result) memcpy(result, a, CALC_BUFFER_SIZE);
- return result;
- }
- if (_desc(b).class == NAN) {
- if (b != result) memcpy(result, b, CALC_BUFFER_SIZE);
- return result;
- }
+ handle_NAN(a, b, result);
- temp = alloca(VALUE_SIZE);
- dividend = alloca(VALUE_SIZE);
+ temp = alloca(value_size);
+ dividend = alloca(value_size);
if (result != a && result != b)
memcpy(&_desc(result), &_desc(a), sizeof(descriptor_t));
_sign(result) = _sign(a) ^ _sign(b);
/* produce nan on 0/0 and inf/inf */
- if (_desc(a).class == ZERO) {
- if (_desc(b).class == ZERO)
+ if (_desc(a).clss == ZERO) {
+ if (_desc(b).clss == ZERO)
/* 0/0 -> nan */
fc_get_qnan(_desc(a).exponent_size, _desc(a).mantissa_size, result);
else
/* 0/x -> a */
- if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, CALC_BUFFER_SIZE-1);
+ if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, calc_buffer_size-1);
return result;
}
- if (_desc(b).class == INF) {
- if (_desc(a).class == INF)
+ if (_desc(b).clss == INF) {
+ if (_desc(a).clss == INF)
/* inf/inf -> nan */
fc_get_qnan(_desc(a).exponent_size, _desc(a).mantissa_size, result);
else {
sc_val_from_ulong(0, NULL);
_save_result(_exp(result));
_save_result(_mant(result));
- _desc(result).class = ZERO;
+ _desc(result).clss = ZERO;
}
return result;
}
- if (_desc(a).class == INF) {
+ if (_desc(a).clss == INF) {
/* inf/x -> inf */
- if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, CALC_BUFFER_SIZE-1);
+ if (a != result) memcpy(result+SIGN_POS+1, a+SIGN_POS+1, calc_buffer_size-1);
return result;
}
- if (_desc(b).class == ZERO) {
+ if (_desc(b).clss == ZERO) {
/* division by zero */
if (_sign(result))
fc_get_minusinf(_desc(a).exponent_size, _desc(a).mantissa_size, result);
sc_add(_exp(result), temp, _exp(result));
/* mixed normal, subnormal values introduce an error of 1, correct it */
- if ((_desc(a).class == SUBNORMAL) ^ (_desc(b).class == SUBNORMAL))
+ if ((_desc(a).clss == SUBNORMAL) ^ (_desc(b).clss == SUBNORMAL))
{
sc_val_from_ulong(1, temp);
sc_add(_exp(result), temp, _exp(result));
_shift_left(_mant(a), temp, dividend);
{
- char *divisor = alloca(CALC_BUFFER_SIZE);
+ char *divisor = alloca(calc_buffer_size);
sc_val_from_ulong(1, divisor);
_shift_right(_mant(b), divisor, divisor);
sc_div(dividend, divisor, _mant(result));
return _normalize(result, result, sc_had_carry());
}
-void _power_of_ten(int exp, descriptor_t *desc, char *result)
+#if 0
+static void _power_of_ten(int exp, descriptor_t *desc, char *result)
{
char *build;
char *temp;
if (desc != NULL)
memcpy(&_desc(result), desc, sizeof(descriptor_t));
- build = alloca(VALUE_SIZE);
- temp = alloca(VALUE_SIZE);
+ build = alloca(value_size);
+ temp = alloca(value_size);
sc_val_from_ulong((1 << _desc(result).exponent_size)/2-1, _exp(result));
}
_save_result(build);
- /* temp is amount of leftshift needed to put the value left of the radix point */
+ /* temp is amount of left shift needed to put the value left of the radix point */
sc_val_from_ulong(_desc(result).mantissa_size + 2, temp);
_shift_left(build, temp, _mant(result));
_normalize(result, result, 0);
}
}
+#endif
+/**
+ * Truncate the fractional part away.
+ *
+ * This does not clip to any integer rang.
+ */
static char* _trunc(const char *a, char *result)
{
- /* when exponent == 0 all bits left of the radix point
+ /*
+ * When exponent == 0 all bits left of the radix point
* are the integral part of the value. For 15bit exp_size
- * this would require a leftshift of max. 16383 bits which
+ * this would require a left shift of max. 16383 bits which
* is too much.
* But it is enough to ensure that no bit right of the radix
* point remains set. This restricts the interesting
* exponents to the interval [0, mant_size-1].
* Outside this interval the truncated value is either 0 or
- * it is are already truncated */
+ * it does not have fractional parts.
+ */
int exp_bias, exp_val;
char *temp;
- temp = alloca(VALUE_SIZE);
+ temp = alloca(value_size);
if (a != result)
memcpy(&_desc(result), &_desc(a), sizeof(descriptor_t));
exp_bias = (1<<_desc(a).exponent_size)/2-1;
- exp_val = sc_val_to_long(_exp(a)) - exp_bias;
+ exp_val = sc_val_to_long(_exp(a)) - exp_bias;
if (exp_val < 0) {
sc_val_from_ulong(0, NULL);
_save_result(_exp(result));
_save_result(_mant(result));
- _desc(result).class = ZERO;
+ _desc(result).clss = ZERO;
return result;
}
if (exp_val > _desc(a).mantissa_size) {
if (a != result)
- memcpy(result, a, CALC_BUFFER_SIZE);
+ memcpy(result, a, calc_buffer_size);
return result;
}
/* and the mask and return the result */
sc_and(_mant(a), temp, _mant(result));
- if (a != result) memcpy(_exp(result), _exp(a), VALUE_SIZE);
+ if (a != result) memcpy(_exp(result), _exp(a), value_size);
return result;
}
/*
* This does value sanity checking(or should do it), sets up any prerequisites,
* calls the proper internal functions, clears up and returns
- * the result */
+ * the result
+ */
char* _calc(const char *a, const char *b, int opcode, char *result)
{
char *temp;
/* make the value with the bigger exponent the first one */
TRACEPRINTF(("+ %s ", fc_print(b, buffer, 100, FC_PACKED)));
if (sc_comp(_exp(a), _exp(b)) == -1)
- _add(b, a, result);
+ _fadd(b, a, result);
else
- _add(a, b, result);
+ _fadd(a, b, result);
break;
case FC_sub:
TRACEPRINTF(("- %s ", fc_print(b, buffer, 100, FC_PACKED)));
- temp = alloca(CALC_BUFFER_SIZE);
- memcpy(temp, b, CALC_BUFFER_SIZE);
+ temp = alloca(calc_buffer_size);
+ memcpy(temp, b, calc_buffer_size);
_sign(temp) = !_sign(b);
if (sc_comp(_exp(a), _exp(temp)) == -1)
- _add(temp, a, result);
+ _fadd(temp, a, result);
else
- _add(a, temp, result);
+ _fadd(a, temp, result);
break;
case FC_mul:
TRACEPRINTF(("* %s ", fc_print(b, buffer, 100, FC_PACKED)));
- _mul(a, b, result);
+ _fmul(a, b, result);
break;
case FC_div:
TRACEPRINTF(("/ %s ", fc_print(b, buffer, 100, FC_PACKED)));
- _div(a, b, result);
+ _fdiv(a, b, result);
break;
case FC_neg:
TRACEPRINTF(("negated "));
- if (a != result) memcpy(result, a, CALC_BUFFER_SIZE);
+ if (a != result) memcpy(result, a, calc_buffer_size);
_sign(result) = !_sign(a);
break;
case FC_int:
+ TRACEPRINTF(("truncated to integer "));
_trunc(a, result);
break;
case FC_rnd:
+ TRACEPRINTF(("rounded to integer "));
+ assert(0);
break;
}
return calc_buffer;
}
-const int fc_get_buffer_length(void)
+int fc_get_buffer_length(void)
{
- return CALC_BUFFER_SIZE;
+ return calc_buffer_size;
}
char* fc_val_from_str(const char *str, unsigned int len, char exp_size, char mant_size, char *result)
if (result == NULL) result = calc_buffer;
- exp_val = alloca(VALUE_SIZE);
- power_val = alloca(CALC_BUFFER_SIZE);
+ exp_val = alloca(value_size);
+ power_val = alloca(calc_buffer_size);
mant_str = alloca((len)?(len):(strlen(str)));
_desc(result).exponent_size = exp_size;
_desc(result).mantissa_size = mant_size;
- _desc(result).class = NORMAL;
+ _desc(result).clss = NORMAL;
old_str = str;
pos = 0;
_power_of_ten(exp_int, &_desc(result), power_val);
- _div(result, power_val, result);
+ _fdiv(result, power_val, result);
return result;
#else
#endif
if (result == NULL) result = calc_buffer;
- temp = alloca(VALUE_SIZE);
+ temp = alloca(value_size);
_desc(result).exponent_size = exp_size;
_desc(result).mantissa_size = mant_size;
/* sign and flag suffice to identify nan or inf, no exponent/mantissa
* encoding is needed. the function can return immediately in these cases */
if (isnan(l)) {
- _desc(result).class = NAN;
+ _desc(result).clss = NAN;
TRACEPRINTF(("val_from_float resulted in NAN\n"));
return result;
}
else if (isinf(l)) {
- _desc(result).class = INF;
+ _desc(result).clss = INF;
TRACEPRINTF(("val_from_float resulted in %sINF\n", (_sign(result)==1)?"-":""));
return result;
}
_normalize(result, result, 0);
- TRACEPRINTF(("val_from_float results in %s\n", fc_print(result, temp, CALC_BUFFER_SIZE, FC_PACKED)));
+ TRACEPRINTF(("val_from_float results in %s\n", fc_print(result, temp, calc_buffer_size, FC_PACKED)));
return result;
}
char result_mantissa = 52;
#endif
- temp = alloca(CALC_BUFFER_SIZE);
+ temp = alloca(calc_buffer_size);
#ifdef HAVE_EXPLICIT_ONE
value = fc_cast(val, result_exponent, result_mantissa-1, temp);
#else
mantissa0 |= sc_sub_bits(_mant(value), result_mantissa, byte_offset) << ((byte_offset-4)<<3);
#ifdef HAVE_LONG_DOUBLE
- mantissa0 &= 0x000FFFFF; /* get rid of garbage */
buildval.val.high = sign << 15;
buildval.val.high |= exponent;
buildval.val.mid = mantissa0;
buildval.val.low = mantissa1;
#else /* no long double */
+ mantissa0 &= 0x000FFFFF; /* get rid of garbage */
buildval.val.high = sign << 31;
buildval.val.high |= exponent << 20;
buildval.val.high |= mantissa0;
int exp_offset, val_bias, res_bias;
if (result == NULL) result = calc_buffer;
- temp = alloca(VALUE_SIZE);
+ temp = alloca(value_size);
if (_desc(value).exponent_size == exp_size && _desc(value).mantissa_size == mant_size)
{
- if (value != result) memcpy(result, value, CALC_BUFFER_SIZE);
+ if (value != result) memcpy(result, value, calc_buffer_size);
return result;
}
/* set the descriptor of the new value */
_desc(result).exponent_size = exp_size;
_desc(result).mantissa_size = mant_size;
- _desc(result).class = _desc(value).class;
+ _desc(result).clss = _desc(value).clss;
_sign(result) = _sign(value);
sc_add(_exp(value), temp, _exp(result));
/* _normalize expects normalized radix point */
- if (_desc(val).class == SUBNORMAL) {
+ if (_desc(val).clss == SUBNORMAL) {
sc_val_from_ulong(1, NULL);
_shift_left(_mant(val), sc_get_buffer(), _mant(result));
} else if (value != result) {
- memcpy(_mant(result), _mant(value), VALUE_SIZE);
+ memcpy(_mant(result), _mant(value), value_size);
} else {
- memmove(_mant(result), _mant(value), VALUE_SIZE);
+ memmove(_mant(result), _mant(value), value_size);
}
_normalize(result, result, 0);
- TRACEPRINTF(("Cast results in %s\n", fc_print(result, temp, VALUE_SIZE, FC_PACKED)));
+ TRACEPRINTF(("Cast results in %s\n", fc_print(result, temp, value_size, FC_PACKED)));
return result;
}
_desc(result).exponent_size = exponent_size;
_desc(result).mantissa_size = mantissa_size;
- _desc(result).class = NORMAL;
+ _desc(result).clss = NORMAL;
_sign(result) = 0;
_desc(result).exponent_size = exponent_size;
_desc(result).mantissa_size = mantissa_size;
- _desc(result).class = NAN;
+ _desc(result).clss = NAN;
_sign(result) = 0;
sc_val_from_ulong((1<<exponent_size)-1, _exp(result));
- /* signalling nan has non-zero mantissa with msb not set */
+ /* signaling NaN has non-zero mantissa with msb not set */
sc_val_from_ulong(1, _mant(result));
return result;
_desc(result).exponent_size = exponent_size;
_desc(result).mantissa_size = mantissa_size;
- _desc(result).class = NAN;
+ _desc(result).clss = NAN;
_sign(result) = 0;
sc_val_from_ulong((1<<exponent_size)-1, _exp(result));
- /* quiet nan has the msb of the mantissa set, so shift one there */
+ /* quiet NaN has the msb of the mantissa set, so shift one there */
sc_val_from_ulong(1, _mant(result));
/* mantissa_size >+< 1 because of two extra rounding bits */
sc_val_from_ulong(mantissa_size + 1, NULL);
_desc(result).exponent_size = exponent_size;
_desc(result).mantissa_size = mantissa_size;
- _desc(result).class = NORMAL;
+ _desc(result).clss = NORMAL;
_sign(result) = 0;
{
const char *val_a = (const char*)a;
const char *val_b = (const char*)b;
+ int mul = 1;
+
+ /*
+ * shortcut: if both values are identical, they are either
+ * Unordered if NaN or equal
+ */
+ if (a == b)
+ return _desc(val_a).clss == NAN ? 2 : 0;
+
+ /* unordered if one is a NaN */
+ if (_desc(val_a).clss == NAN || _desc(val_b).clss == NAN)
+ return 2;
- /* unordered */
- if (_desc(val_a).class == NAN || _desc(val_b).class == NAN) return 2;
/* zero is equal independent of sign */
- if ((_desc(val_a).class == ZERO) && (_desc(val_b).class == ZERO)) return 0;
+ if ((_desc(val_a).clss == ZERO) && (_desc(val_b).clss == ZERO))
+ return 0;
+
/* different signs make compare easy */
- if (_sign(val_a) != _sign(val_b)) return (_sign(val_a)==0)?(1):(-1);
+ if (_sign(val_a) != _sign(val_b))
+ return (_sign(val_a)==0)?(1):(-1);
+
+ mul = _sign(a) ? -1 : 1;
+
/* both infinity means equality */
- if ((_desc(val_a).class == INF) && (_desc(val_b).class == INF)) return 0;
+ if ((_desc(val_a).clss == INF) && (_desc(val_b).clss == INF))
+ return 0;
+
/* infinity is bigger than the rest */
- if (_desc(val_a).class == INF) return _sign(val_a)?(-1):(1);
- if (_desc(val_b).class == INF) return _sign(val_b)?(1):(-1);
+ if (_desc(val_a).clss == INF)
+ return 1 * mul;
+ if (_desc(val_b).clss == INF)
+ return -1 * mul;
+ /* check first exponent, that mantissa if equal */
switch (sc_comp(_exp(val_a), _exp(val_b))) {
case -1:
- return -1;
+ return -1 * mul;
case 1:
- return 1;
+ return 1 * mul;
case 0:
- return sc_comp(_mant(val_a), _mant(val_b));
+ return sc_comp(_mant(val_a), _mant(val_b)) * mul;
default:
return 2;
}
int fc_is_zero(const void *a)
{
- return _desc((const char*)a).class == ZERO;
+ return _desc(a).clss == ZERO;
}
int fc_is_negative(const void *a)
{
- return _sign((const char*)a);
+ return _sign(a);
}
int fc_is_inf(const void *a)
{
- return _desc(a).class == INF;
+ return _desc(a).clss == INF;
}
int fc_is_nan(const void *a)
{
- return _desc(a).class == NAN;
+ return _desc(a).clss == NAN;
}
int fc_is_subnormal(const void *a)
{
- return _desc(a).class == SUBNORMAL;
+ return _desc(a).clss == SUBNORMAL;
}
char *fc_print(const void *a, char *buf, int buflen, unsigned base)
val = (const char*)a;
- mul_1 = alloca(CALC_BUFFER_SIZE);
+ mul_1 = alloca(calc_buffer_size);
switch (base) {
case FC_DEC:
- switch (_desc(val).class) {
+ switch (_desc(val).clss) {
case INF:
if (buflen >= 8+_sign(val)) sprintf(buf, "%sINFINITY", _sign(val)?"-":"");
else snprintf(buf, buflen, "%sINF", _sign(val)?"-":NULL);
break;
case FC_HEX:
- switch (_desc(val).class) {
+ switch (_desc(val).clss) {
case INF:
if (buflen >= 8+_sign(val)) sprintf(buf, "%sINFINITY", _sign(val)?"-":"");
else snprintf(buf, buflen, "%sINF", _sign(val)?"-":NULL);
case FC_PACKED:
default:
- snprintf(buf, buflen, "%s", sc_print(_pack(val, mul_1), VALUE_SIZE*4, SC_HEX));
+ snprintf(buf, buflen, "%s", sc_print(_pack(val, mul_1), value_size*4, SC_HEX, 0));
+ buf[buflen - 1] = '\0';
break;
}
return buf;
/* this is used to cache the packed version of the value */
static char *pack = NULL;
- if (pack == NULL) pack = malloc(VALUE_SIZE);
+ if (pack == NULL) pack = xmalloc(value_size);
if (value != NULL)
_pack((const char*)value, pack);
fc_rounding_mode_t fc_set_rounding_mode(fc_rounding_mode_t mode)
{
if (mode == FC_TONEAREST || mode == FC_TOPOSITIVE || mode == FC_TONEGATIVE || mode == FC_TOZERO)
- ROUNDING_MODE = mode;
+ rounding_mode = mode;
- return ROUNDING_MODE;
+ return rounding_mode;
}
fc_rounding_mode_t fc_get_rounding_mode(void)
{
- return ROUNDING_MODE;
+ return rounding_mode;
}
void init_fltcalc(int precision)
if (max_precision < precision)
printf("WARING: not enough precision available, using %d\n", max_precision);
- ROUNDING_MODE = FC_TONEAREST;
- VALUE_SIZE = sc_get_buffer_length();
- SIGN_POS = 0;
- EXPONENT_POS = SIGN_POS + sizeof(char);
- MANTISSA_POS = EXPONENT_POS + VALUE_SIZE;
- DESCRIPTOR_POS = MANTISSA_POS + VALUE_SIZE;
- CALC_BUFFER_SIZE = DESCRIPTOR_POS + sizeof(descriptor_t);
-
- calc_buffer = malloc(CALC_BUFFER_SIZE);
- DEBUGPRINTF(("init fltcalc:\n\tVALUE_SIZE = %d\n\tSIGN_POS = %d\n\tEXPONENT_POS = %d\n\tMANTISSA_POS = %d\n\tDESCRIPTOR_POS = %d\n\tCALC_BUFFER_SIZE = %d\n\tcalc_buffer = %p\n\n", VALUE_SIZE, SIGN_POS, EXPONENT_POS, MANTISSA_POS, DESCRIPTOR_POS, CALC_BUFFER_SIZE, calc_buffer));
+ rounding_mode = FC_TONEAREST;
+ value_size = sc_get_buffer_length();
+ SIGN_POS = 0;
+ EXPONENT_POS = SIGN_POS + sizeof(char);
+ MANTISSA_POS = EXPONENT_POS + value_size;
+ DESCRIPTOR_POS = MANTISSA_POS + value_size;
+ calc_buffer_size = DESCRIPTOR_POS + sizeof(descriptor_t);
+
+ calc_buffer = xmalloc(calc_buffer_size);
+ memset(calc_buffer, 0, calc_buffer_size);
+ DEBUGPRINTF(("init fltcalc:\n\tVALUE_SIZE = %d\n\tSIGN_POS = %d\n\tEXPONENT_POS = %d\n\tMANTISSA_POS = %d\n\tDESCRIPTOR_POS = %d\n\tCALC_BUFFER_SIZE = %d\n\tcalc_buffer = %p\n\n", value_size, SIGN_POS, EXPONENT_POS, MANTISSA_POS, DESCRIPTOR_POS, calc_buffer_size, calc_buffer));
#ifdef HAVE_LONG_DOUBLE
DEBUGPRINTF(("\tUsing long double (1-15-64) interface\n"));
#else