X-Git-Url: http://nsz.repo.hu/git/?a=blobdiff_plain;f=ir%2Ftv%2Ffltcalc.c;h=7c5d6fdd15a1eed8f34c49b0d1db5674254e7965;hb=e70c22e8aa4cc93e7d72f58453c03cb9bd797e31;hp=bc9d803bd1bf9d2a13b7fe372cf0e879e8828d81;hpb=c569340a123b51f08b8e95a2093f4f8aa1bcc27d;p=libfirm diff --git a/ir/tv/fltcalc.c b/ir/tv/fltcalc.c index bc9d803bd..7c5d6fdd1 100644 --- a/ir/tv/fltcalc.c +++ b/ir/tv/fltcalc.c @@ -26,19 +26,24 @@ #endif #ifdef HAVE_INTTYPES_H -#include +# include +#endif +#ifdef HAVE_STRING_H +# include +#endif +#ifdef HAVE_STDLIB_H +# include #endif - -#include -#include -#include -#include #ifdef HAVE_ALLOCA_H # include #endif #ifdef HAVE_MALLOC_H # include #endif +#include +#include + +#include "xmalloc.h" typedef uint32_t UINT32; @@ -82,54 +87,54 @@ typedef union { #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 @@ -147,10 +152,10 @@ typedef struct { 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; @@ -174,25 +179,25 @@ static void _fail_char(const char *str, unsigned int len, int 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; @@ -235,7 +240,7 @@ char* _normalize(const char *in_val, char *out_val, int sticky) 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; @@ -246,9 +251,9 @@ char* _normalize(const char *in_val, char *out_val, int sticky) 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)); @@ -291,7 +296,7 @@ char* _normalize(const char *in_val, char *out_val, int sticky) /* 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 @@ -301,7 +306,7 @@ char* _normalize(const char *in_val, char *out_val, int sticky) 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 @@ -340,26 +345,26 @@ char* _normalize(const char *in_val, char *out_val, int sticky) } /* 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 @@ -389,10 +394,10 @@ char* _normalize(const char *in_val, char *out_val, int sticky) 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: @@ -401,10 +406,10 @@ char* _normalize(const char *in_val, char *out_val, int sticky) } } 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: @@ -417,10 +422,26 @@ char* _normalize(const char *in_val, char *out_val, int sticky) 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; @@ -428,14 +449,7 @@ static char* _add(const char* a, const char* b, char* result) 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) @@ -444,12 +458,12 @@ static char* _add(const char* a, const char* b, char* result) /* 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); @@ -461,50 +475,45 @@ static char* _add(const char* a, const char* b, char* result) 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); @@ -518,7 +527,7 @@ static char* _add(const char* a, const char* b, char* result) /* 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); } @@ -534,60 +543,56 @@ static char* _add(const char* a, const char* b, char* result) /* _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; } @@ -598,7 +603,7 @@ static char* _mul(const char* a, const char* b, char* 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)); @@ -618,21 +623,17 @@ static char* _mul(const char* a, const char* b, char* 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)); @@ -640,18 +641,18 @@ static char* _div(const char* a, const char* b, char* result) _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 { @@ -659,17 +660,17 @@ static char* _div(const char* a, const char* b, char* result) 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); @@ -684,7 +685,7 @@ static char* _div(const char* a, const char* b, char* 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)); @@ -701,7 +702,7 @@ static char* _div(const char* a, const char* b, char* 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)); @@ -710,7 +711,7 @@ static char* _div(const char* a, const char* b, char* result) return _normalize(result, result, sc_had_carry()); } -void _power_of_ten(int exp, descriptor_t *desc, char *result) +static void _power_of_ten(int exp, descriptor_t *desc, char *result) { char *build; char *temp; @@ -722,8 +723,8 @@ void _power_of_ten(int exp, descriptor_t *desc, char *result) 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)); @@ -739,7 +740,7 @@ void _power_of_ten(int exp, descriptor_t *desc, char *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)); @@ -748,41 +749,48 @@ void _power_of_ten(int exp, descriptor_t *desc, char *result) } } +/** + * 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; } @@ -796,7 +804,7 @@ static char* _trunc(const char *a, char *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; } @@ -804,7 +812,8 @@ static char* _trunc(const char *a, char *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; @@ -823,37 +832,40 @@ char* _calc(const char *a, const char *b, int opcode, char *result) /* 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; } @@ -869,9 +881,9 @@ const void *fc_get_buffer(void) 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) @@ -896,13 +908,13 @@ char* fc_val_from_str(const char *str, unsigned int len, char exp_size, char man 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; @@ -1051,7 +1063,7 @@ done: _power_of_ten(exp_int, &_desc(result), power_val); - _div(result, power_val, result); + _fdiv(result, power_val, result); return result; #else @@ -1104,7 +1116,7 @@ char* fc_val_from_float(LLDBL l, char exp_size, char mant_size, char* result) #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; @@ -1115,12 +1127,12 @@ char* fc_val_from_float(LLDBL l, char exp_size, char mant_size, char* result) /* 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; } @@ -1174,7 +1186,7 @@ char* fc_val_from_float(LLDBL l, char exp_size, char mant_size, char* 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; } @@ -1201,7 +1213,7 @@ LLDBL fc_val_to_float(const void *val) 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 @@ -1226,16 +1238,13 @@ LLDBL fc_val_to_float(const void *val) for (; (byte_offset<<3) < result_mantissa; byte_offset++) mantissa0 |= sc_sub_bits(_mant(value), result_mantissa, byte_offset) << ((byte_offset-4)<<3); -#ifndef HAVE_LONG_DOUBLE - mantissa0 &= 0x000FFFFF; /* get rid of garbage */ -#endif - #ifdef HAVE_LONG_DOUBLE 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; @@ -1253,18 +1262,18 @@ char* fc_cast(const void *val, char exp_size, char mant_size, char *result) 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); @@ -1279,17 +1288,17 @@ char* fc_cast(const void *val, char exp_size, char mant_size, char *result) 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; } @@ -1299,7 +1308,7 @@ char* fc_get_max(unsigned int exponent_size, unsigned int mantissa_size, char* r _desc(result).exponent_size = exponent_size; _desc(result).mantissa_size = mantissa_size; - _desc(result).class = NORMAL; + _desc(result).clss = NORMAL; _sign(result) = 0; @@ -1328,13 +1337,13 @@ char* fc_get_snan(unsigned int exponent_size, unsigned int mantissa_size, char * _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<+< 1 because of two extra rounding bits */ sc_val_from_ulong(mantissa_size + 1, NULL); @@ -1367,7 +1376,7 @@ char* fc_get_plusinf(unsigned int exponent_size, unsigned int mantissa_size, cha _desc(result).exponent_size = exponent_size; _desc(result).mantissa_size = mantissa_size; - _desc(result).class = NORMAL; + _desc(result).clss = NORMAL; _sign(result) = 0; @@ -1392,26 +1401,47 @@ int fc_comp(const void *a, const void *b) { 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; } @@ -1419,27 +1449,27 @@ int fc_comp(const void *a, const void *b) 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) @@ -1449,11 +1479,11 @@ 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); @@ -1476,7 +1506,7 @@ char *fc_print(const void *a, char *buf, int buflen, unsigned base) 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); @@ -1498,7 +1528,7 @@ char *fc_print(const void *a, char *buf, int buflen, unsigned base) 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)); break; } return buf; @@ -1509,7 +1539,7 @@ unsigned char fc_sub_bits(const void *value, unsigned num_bits, unsigned byte_of /* 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); @@ -1520,14 +1550,14 @@ unsigned char fc_sub_bits(const void *value, unsigned num_bits, unsigned byte_of 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) @@ -1544,16 +1574,16 @@ 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); + 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)); + calc_buffer = xmalloc(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