Update of /cvsroot/bochs/bochs/fpu
In directory sc8-pr-cvs1.sourceforge.net:/tmp/cvs-serv11688/fpu
Modified Files:
Tag: softfloat-fpu-implementation_ver4_branch
Makefile.in div_Xsig.c softfloat-macros.h softfloat.h
Added Files:
Tag: softfloat-fpu-implementation_ver4_branch
poly_eval.cc softfloat128.cc
Log Message:
Add float128 which will be used later in polynomial approximations evaluation
--- NEW FILE: poly_eval.cc ---
/*============================================================================
* Written for Bochs (x86 achitecture simulator) by
* Stanislav Shwartsman (gate at fidonet.org.il)
* ==========================================================================*/
#define FLOAT128
#include "softfloat.h"
/* arithmetical polynomial evauation */
float128 poly_eval(float128 x, float128 *carr, int n, float_status_t &status)
{
float128 result = carr[n];
for(int i=n-1; i >= 0; i--)
{
result = float128_mul(result, x, status);
result = float128_add(result, carr[i], status);
}
return result;
}
--- NEW FILE: softfloat128.cc ---
/*============================================================================
This C source file is part of the SoftFloat IEC/IEEE Floating-point Arithmetic
Package, Release 2b.
Written by John R. Hauser. This work was made possible in part by the
International Computer Science Institute, located at Suite 600, 1947 Center
Street, Berkeley, California 94704. Funding was partially provided by the
National Science Foundation under grant MIP-9311980. The original version
of this code was written as part of a project to build a fixed-point vector
processor in collaboration with the University of California at Berkeley,
overseen by Profs. Nelson Morgan and John Wawrzynek. More information
is available through the Web page `http://www.cs.berkeley.edu/~jhauser/
arithmetic/SoftFloat.html'.
THIS SOFTWARE IS DISTRIBUTED AS IS, FOR FREE. Although reasonable effort has
been made to avoid it, THIS SOFTWARE MAY CONTAIN FAULTS THAT WILL AT TIMES
RESULT IN INCORRECT BEHAVIOR. USE OF THIS SOFTWARE IS RESTRICTED TO PERSONS
AND ORGANIZATIONS WHO CAN AND WILL TAKE FULL RESPONSIBILITY FOR ALL LOSSES,
COSTS, OR OTHER PROBLEMS THEY INCUR DUE TO THE SOFTWARE, AND WHO FURTHERMORE
EFFECTIVELY INDEMNIFY JOHN HAUSER AND THE INTERNATIONAL COMPUTER SCIENCE
INSTITUTE (possibly via similar legal warning) AGAINST ALL LOSSES, COSTS, OR
OTHER PROBLEMS INCURRED BY THEIR CUSTOMERS AND CLIENTS DUE TO THE SOFTWARE.
Derivative works are acceptable, even for commercial purposes, so long as
(1) the source code for the derivative work includes prominent notice that
the work is derivative, and (2) the source code includes prominent notice with
these four paragraphs for those parts of this code that are retained.
=============================================================================*/
/*============================================================================
* Adapted for Bochs (x86 achitecture simulator) by
* Stanislav Shwartsman (gate at fidonet.org.il)
* ==========================================================================*/
#define FLOAT128
#include "softfloat.h"
#include "softfloat-macros.h"
#include "softfloat-round-pack.h"
#include "softfloat-specialize.h"
#ifdef FLOAT128
/*----------------------------------------------------------------------------
| The pattern for a default generated quadruple-precision NaN. The `high' and
| `low' values hold the most- and least-significant bits, respectively.
*----------------------------------------------------------------------------*/
#define float128_default_nan_hi BX_CONST64(0xFFFF800000000000)
#define float128_default_nan_lo BX_CONST64(0x0000000000000000)
#ifdef BX_BIG_ENDIAN
static const float128 float128_default_nan = { float128_default_nan_hi, float128_default_nan_lo };
#else
static const float128 float128_default_nan = { float128_default_nan_lo, float128_default_nan_hi };
#endif
/*----------------------------------------------------------------------------
| Returns the least-significant 64 fraction bits of the quadruple-precision
| floating-point value `a'.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE Bit64u extractFloat128Frac1(float128 a)
{
return a.lo;
}
/*----------------------------------------------------------------------------
| Returns the most-significant 48 fraction bits of the quadruple-precision
| floating-point value `a'.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE Bit64u extractFloat128Frac0(float128 a)
{
return a.hi & BX_CONST64(0x0000FFFFFFFFFFFF);
}
/*----------------------------------------------------------------------------
| Returns the exponent bits of the quadruple-precision floating-point value
| `a'.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE Bit32s extractFloat128Exp(float128 a)
{
return (a.hi>>48) & 0x7FFF;
}
/*----------------------------------------------------------------------------
| Returns the sign bit of the quadruple-precision floating-point value `a'.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE int extractFloat128Sign(float128 a)
{
return a.hi >> 63;
}
/*----------------------------------------------------------------------------
| Returns 1 if the quadruple-precision floating-point value `a' is a NaN;
| otherwise returns 0.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE int float128_is_nan(float128 a)
{
return (BX_CONST64(0xFFFE000000000000) <= (Bit64u) (a.hi<<1))
&& (a.lo || (a.hi & BX_CONST64(0x0000FFFFFFFFFFFF)));
}
/*----------------------------------------------------------------------------
| Returns 1 if the quadruple-precision floating-point value `a' is a
| signaling NaN; otherwise returns 0.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE int float128_is_signaling_nan(float128 a)
{
return (((a.hi>>47) & 0xFFFF) == 0xFFFE)
&& (a.lo || (a.hi & BX_CONST64(0x00007FFFFFFFFFFF)));
}
/*----------------------------------------------------------------------------
| Returns the result of converting the quadruple-precision floating-point NaN
| `a' to the canonical NaN format. If `a' is a signaling NaN, the invalid
| exception is raised.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE commonNaNT float128ToCommonNaN(float128 a, float_status_t &status)
{
commonNaNT z;
if (float128_is_signaling_nan(a)) float_raise(status, float_flag_invalid);
z.sign = a.hi>>63;
shortShift128Left(a.hi, a.lo, 16, &z.hi, &z.lo);
return z;
}
/*----------------------------------------------------------------------------
| Returns the result of converting the canonical NaN `a' to the quadruple-
| precision floating-point format.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE float128 commonNaNToFloat128(commonNaNT a)
{
float128 z;
shift128Right(a.hi, a.lo, 16, &z.hi, &z.lo);
z.hi |= (((Bit64u) a.sign)<<63) | BX_CONST64(0x7FFF800000000000);
return z;
}
/*----------------------------------------------------------------------------
| Takes two quadruple-precision floating-point values `a' and `b', one of
| which is a NaN, and returns the appropriate NaN result. If either `a' or
| `b' is a signaling NaN, the invalid exception is raised.
*----------------------------------------------------------------------------*/
static float128 propagateFloat128NaN(float128 a, float128 b, float_status_t &status)
{
int aIsNaN, aIsSignalingNaN, bIsNaN, bIsSignalingNaN;
aIsNaN = float128_is_nan(a);
aIsSignalingNaN = float128_is_signaling_nan(a);
bIsNaN = float128_is_nan(b);
bIsSignalingNaN = float128_is_signaling_nan(b);
a.hi |= BX_CONST64( 0x0000800000000000 );
b.hi |= BX_CONST64( 0x0000800000000000 );
if (aIsSignalingNaN | bIsSignalingNaN) float_raise(status, float_flag_invalid);
if (aIsSignalingNaN) {
if (bIsSignalingNaN) goto returnLargerSignificand;
return bIsNaN ? b : a;
}
else if (aIsNaN) {
if (bIsSignalingNaN | !bIsNaN) return a;
returnLargerSignificand:
if (lt128(a.hi<<1, a.lo, b.hi<<1, b.lo)) return b;
if (lt128(b.hi<<1, b.lo, a.hi<<1, a.lo)) return a;
return (a.hi < b.hi) ? a : b;
}
else {
return b;
}
}
/*----------------------------------------------------------------------------
| Normalizes the subnormal quadruple-precision floating-point value
| represented by the denormalized significand formed by the concatenation of
| `aSig0' and `aSig1'. The normalized exponent is stored at the location
| pointed to by `zExpPtr'. The most significant 49 bits of the normalized
| significand are stored at the location pointed to by `zSig0Ptr', and the
| least significant 64 bits of the normalized significand are stored at the
| location pointed to by `zSig1Ptr'.
*----------------------------------------------------------------------------*/
static void
normalizeFloat128Subnormal(
Bit64u aSig0,
Bit64u aSig1,
Bit32s *zExpPtr,
Bit64u *zSig0Ptr,
Bit64u *zSig1Ptr
)
{
int shiftCount;
if (aSig0 == 0) {
shiftCount = countLeadingZeros64(aSig1) - 15;
if (shiftCount < 0) {
*zSig0Ptr = aSig1 >>(-shiftCount);
*zSig1Ptr = aSig1 << (shiftCount & 63);
}
else {
*zSig0Ptr = aSig1 << shiftCount;
*zSig1Ptr = 0;
}
*zExpPtr = - shiftCount - 63;
}
else {
shiftCount = countLeadingZeros64(aSig0) - 15;
shortShift128Left(aSig0, aSig1, shiftCount, zSig0Ptr, zSig1Ptr);
*zExpPtr = 1 - shiftCount;
}
}
/*----------------------------------------------------------------------------
| Packs the sign `zSign', the exponent `zExp', and the significand formed
| by the concatenation of `zSig0' and `zSig1' into a quadruple-precision
| floating-point value, returning the result. After being shifted into the
| proper positions, the three fields `zSign', `zExp', and `zSig0' are simply
| added together to form the most significant 32 bits of the result. This
| means that any integer portion of `zSig0' will be added into the exponent.
| Since a properly normalized significand will have an integer portion equal
| to 1, the `zExp' input should be 1 less than the desired result exponent
| whenever `zSig0' and `zSig1' concatenated form a complete, normalized
| significand.
*----------------------------------------------------------------------------*/
BX_CPP_INLINE float128 packFloat128(int zSign, Bit32s zExp, Bit64u zSig0, Bit64u zSig1)
{
float128 z;
z.lo = zSig1;
z.hi = (((Bit64u) zSign)<<63) + (((Bit64u) zExp)<<48) + zSig0;
return z;
}
/*----------------------------------------------------------------------------
| Takes an abstract floating-point value having sign `zSign', exponent `zExp',
| and extended significand formed by the concatenation of `zSig0', `zSig1',
| and `zSig2', and returns the proper quadruple-precision floating-point value
| corresponding to the abstract input. Ordinarily, the abstract value is
| simply rounded and packed into the quadruple-precision format, with the
| inexact exception raised if the abstract input cannot be represented
| exactly. However, if the abstract value is too large, the overflow and
| inexact exceptions are raised and an infinity or maximal finite value is
| returned. If the abstract value is too small, the input value is rounded to
| a subnormal number, and the underflow and inexact exceptions are raised if
| the abstract input cannot be represented exactly as a subnormal quadruple-
| precision floating-point number.
| The input significand must be normalized or smaller. If the input
| significand is not normalized, `zExp' must be 0; in that case, the result
| returned is a subnormal number, and it must not require rounding. In the
| usual case that the input significand is normalized, `zExp' must be 1 less
| than the ``true'' floating-point exponent. The handling of underflow and
| overflow follows the IEC/IEEE Standard for Binary Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
static float128
roundAndPackFloat128(
int zSign, Bit32s zExp, Bit64u zSig0, Bit64u zSig1, Bit64u zSig2, float_status_t &status)
{
int increment = ((Bit64s) zSig2 < 0);
if (0x7FFD <= (Bit32u) zExp) {
if ((0x7FFD < zExp)
|| ((zExp == 0x7FFD)
&& eq128(BX_CONST64(0x0001FFFFFFFFFFFF),
BX_CONST64(0xFFFFFFFFFFFFFFFF), zSig0, zSig1)
&& increment))
{
float_raise(status, float_flag_overflow | float_flag_inexact);
return packFloat128(zSign, 0x7FFF, 0, 0);
}
if (zExp < 0) {
int isTiny = (zExp < -1)
|| ! increment
|| lt128(zSig0, zSig1,
BX_CONST64(0x0001FFFFFFFFFFFF),
BX_CONST64(0xFFFFFFFFFFFFFFFF)
);
shift128ExtraRightJamming(
zSig0, zSig1, zSig2, -zExp, &zSig0, &zSig1, &zSig2);
zExp = 0;
if (isTiny && zSig2) float_raise(status, float_flag_underflow);
increment = ((Bit64s) zSig2 < 0);
}
}
if (zSig2) float_raise(status, float_flag_inexact);
if (increment) {
add128(zSig0, zSig1, 0, 1, &zSig0, &zSig1);
zSig1 &= ~((zSig2 + zSig2 == 0) & 1);
}
else {
if ((zSig0 | zSig1) == 0) zExp = 0;
}
return packFloat128(zSign, zExp, zSig0, zSig1);
}
/*----------------------------------------------------------------------------
| Takes an abstract floating-point value having sign `zSign', exponent `zExp',
| and significand formed by the concatenation of `zSig0' and `zSig1', and
| returns the proper quadruple-precision floating-point value corresponding
| to the abstract input. This routine is just like `roundAndPackFloat128'
| except that the input significand has fewer bits and does not have to be
| normalized. In all cases, `zExp' must be 1 less than the ``true'' floating-
| point exponent.
*----------------------------------------------------------------------------*/
static float128
normalizeRoundAndPackFloat128(
int zSign, Bit32s zExp, Bit64u zSig0, Bit64u zSig1, float_status_t &status)
{
Bit64u zSig2;
if (zSig0 == 0) {
zSig0 = zSig1;
zSig1 = 0;
zExp -= 64;
}
int shiftCount = countLeadingZeros64(zSig0) - 15;
if (0 <= shiftCount) {
zSig2 = 0;
shortShift128Left(zSig0, zSig1, shiftCount, &zSig0, &zSig1);
}
else {
shift128ExtraRightJamming(
zSig0, zSig1, 0, -shiftCount, &zSig0, &zSig1, &zSig2);
}
zExp -= shiftCount;
return roundAndPackFloat128(zSign, zExp, zSig0, zSig1, zSig2, status);
}
/*----------------------------------------------------------------------------
| Returns the result of converting the extended double-precision floating-
| point value `a' to the quadruple-precision floating-point format. The
| conversion is performed according to the IEC/IEEE Standard for Binary
| Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
float128 floatx80_to_float128(floatx80 a, float_status_t &status)
{
Bit64u zSig0, zSig1;
Bit64u aSig = extractFloatx80Frac(a);
Bit32s aExp = extractFloatx80Exp(a);
int aSign = extractFloatx80Sign(a);
if ((aExp == 0x7FFF) && (Bit64u) (aSig<<1))
return commonNaNToFloat128(floatx80ToCommonNaN(a, status));
shift128Right(aSig<<1, 0, 16, &zSig0, &zSig1);
return packFloat128(aSign, aExp, zSig0, zSig1);
}
/*----------------------------------------------------------------------------
| Returns the result of converting the quadruple-precision floating-point
| value `a' to the extended double-precision floating-point format. The
| conversion is performed according to the IEC/IEEE Standard for Binary
| Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
floatx80 float128_to_floatx80(float128 a, float_status_t &status)
{
Bit32s aExp;
Bit64u aSig0, aSig1;
aSig1 = extractFloat128Frac1(a);
aSig0 = extractFloat128Frac0(a);
aExp = extractFloat128Exp(a);
int aSign = extractFloat128Sign(a);
if (aExp == 0x7FFF) {
if (aSig0 | aSig1)
return commonNaNToFloatx80(float128ToCommonNaN(a, status));
return packFloatx80(aSign, 0x7FFF, BX_CONST64(0x8000000000000000));
}
if (aExp == 0) {
if ((aSig0 | aSig1) == 0) return packFloatx80(aSign, 0, 0);
normalizeFloat128Subnormal(aSig0, aSig1, &aExp, &aSig0, &aSig1);
}
else aSig0 |= BX_CONST64(0x0001000000000000);
shortShift128Left(aSig0, aSig1, 15, &aSig0, &aSig1);
return roundAndPackFloatx80(80, aSign, aExp, aSig0, aSig1, status);
}
/*----------------------------------------------------------------------------
| Returns the result of adding the absolute values of the quadruple-precision
| floating-point values `a' and `b'. If `zSign' is 1, the sum is negated
| before being returned. `zSign' is ignored if the result is a NaN.
| The addition is performed according to the IEC/IEEE Standard for Binary
| Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
static float128 addFloat128Sigs(float128 a, float128 b, int zSign, float_status_t &status)
{
Bit32s aExp, bExp, zExp;
Bit64u aSig0, aSig1, bSig0, bSig1, zSig0, zSig1, zSig2;
Bit32s expDiff;
aSig1 = extractFloat128Frac1(a);
aSig0 = extractFloat128Frac0(a);
aExp = extractFloat128Exp(a);
bSig1 = extractFloat128Frac1(b);
bSig0 = extractFloat128Frac0(b);
bExp = extractFloat128Exp(b);
expDiff = aExp - bExp;
if (0 < expDiff) {
if (aExp == 0x7FFF) {
if (aSig0 | aSig1) return propagateFloat128NaN(a, b, status);
return a;
}
if (bExp == 0) --expDiff;
else {
bSig0 |= BX_CONST64(0x0001000000000000);
}
shift128ExtraRightJamming(
bSig0, bSig1, 0, expDiff, &bSig0, &bSig1, &zSig2);
zExp = aExp;
}
else if (expDiff < 0) {
if (bExp == 0x7FFF) {
if (bSig0 | bSig1) return propagateFloat128NaN(a, b, status);
return packFloat128(zSign, 0x7FFF, 0, 0);
}
if (aExp == 0) ++expDiff;
else {
aSig0 |= BX_CONST64(0x0001000000000000);
}
shift128ExtraRightJamming(
aSig0, aSig1, 0, - expDiff, &aSig0, &aSig1, &zSig2);
zExp = bExp;
}
else {
if (aExp == 0x7FFF) {
if (aSig0 | aSig1 | bSig0 | bSig1)
return propagateFloat128NaN(a, b, status);
return a;
}
add128(aSig0, aSig1, bSig0, bSig1, &zSig0, &zSig1);
if (aExp == 0) return packFloat128(zSign, 0, zSig0, zSig1);
zSig2 = 0;
zSig0 |= BX_CONST64(0x0002000000000000);
zExp = aExp;
goto shiftRight1;
}
aSig0 |= BX_CONST64(0x0001000000000000);
add128(aSig0, aSig1, bSig0, bSig1, &zSig0, &zSig1);
--zExp;
if (zSig0 < BX_CONST64(0x0002000000000000)) goto roundAndPack;
++zExp;
shiftRight1:
shift128ExtraRightJamming(zSig0, zSig1, zSig2, 1, &zSig0, &zSig1, &zSig2);
roundAndPack:
return roundAndPackFloat128(zSign, zExp, zSig0, zSig1, zSig2, status);
}
/*----------------------------------------------------------------------------
| Returns the result of subtracting the absolute values of the quadruple-
| precision floating-point values `a' and `b'. If `zSign' is 1, the
| difference is negated before being returned. `zSign' is ignored if the
| result is a NaN. The subtraction is performed according to the IEC/IEEE
| Standard for Binary Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
static float128 subFloat128Sigs(float128 a, float128 b, int zSign, float_status_t &status)
{
Bit32s aExp, bExp, zExp;
Bit64u aSig0, aSig1, bSig0, bSig1, zSig0, zSig1;
Bit32s expDiff;
aSig1 = extractFloat128Frac1(a);
aSig0 = extractFloat128Frac0(a);
aExp = extractFloat128Exp(a);
bSig1 = extractFloat128Frac1(b);
bSig0 = extractFloat128Frac0(b);
bExp = extractFloat128Exp(b);
expDiff = aExp - bExp;
shortShift128Left(aSig0, aSig1, 14, &aSig0, &aSig1);
shortShift128Left(bSig0, bSig1, 14, &bSig0, &bSig1);
if (0 < expDiff) goto aExpBigger;
if (expDiff < 0) goto bExpBigger;
if (aExp == 0x7FFF) {
if (aSig0 | aSig1 | bSig0 | bSig1)
return propagateFloat128NaN(a, b, status);
float_raise(status, float_flag_invalid);
return float128_default_nan;
}
if (aExp == 0) {
aExp = 1;
bExp = 1;
}
if (bSig0 < aSig0) goto aBigger;
if (aSig0 < bSig0) goto bBigger;
if (bSig1 < aSig1) goto aBigger;
if (aSig1 < bSig1) goto bBigger;
return packFloat128(0, 0, 0, 0);
bExpBigger:
if (bExp == 0x7FFF) {
if (bSig0 | bSig1) return propagateFloat128NaN(a, b, status);
return packFloat128(zSign ^ 1, 0x7FFF, 0, 0);
}
if (aExp == 0) ++expDiff;
else {
aSig0 |= BX_CONST64(0x4000000000000000);
}
shift128RightJamming(aSig0, aSig1, - expDiff, &aSig0, &aSig1);
bSig0 |= BX_CONST64(0x4000000000000000);
bBigger:
sub128(bSig0, bSig1, aSig0, aSig1, &zSig0, &zSig1);
zExp = bExp;
zSign ^= 1;
goto normalizeRoundAndPack;
aExpBigger:
if (aExp == 0x7FFF) {
if (aSig0 | aSig1) return propagateFloat128NaN(a, b, status);
return a;
}
if (bExp == 0) --expDiff;
else {
bSig0 |= BX_CONST64(0x4000000000000000);
}
shift128RightJamming(bSig0, bSig1, expDiff, &bSig0, &bSig1);
aSig0 |= BX_CONST64(0x4000000000000000);
aBigger:
sub128(aSig0, aSig1, bSig0, bSig1, &zSig0, &zSig1);
zExp = aExp;
normalizeRoundAndPack:
--zExp;
return normalizeRoundAndPackFloat128(zSign, zExp - 14, zSig0, zSig1, status);
}
/*----------------------------------------------------------------------------
| Returns the result of adding the quadruple-precision floating-point values
| `a' and `b'. The operation is performed according to the IEC/IEEE Standard
| for Binary Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
float128 float128_add(float128 a, float128 b, float_status_t &status)
{
int aSign = extractFloat128Sign(a);
int bSign = extractFloat128Sign(b);
if (aSign == bSign) {
return addFloat128Sigs(a, b, aSign, status);
}
else {
return subFloat128Sigs(a, b, aSign, status);
}
}
/*----------------------------------------------------------------------------
| Returns the result of subtracting the quadruple-precision floating-point
| values `a' and `b'. The operation is performed according to the IEC/IEEE
| Standard for Binary Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
float128 float128_sub(float128 a, float128 b, float_status_t &status)
{
int aSign = extractFloat128Sign(a);
int bSign = extractFloat128Sign(b);
if (aSign == bSign) {
return subFloat128Sigs(a, b, aSign, status);
}
else {
return addFloat128Sigs(a, b, aSign, status);
}
}
/*----------------------------------------------------------------------------
| Returns the result of multiplying the quadruple-precision floating-point
| values `a' and `b'. The operation is performed according to the IEC/IEEE
| Standard for Binary Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
float128 float128_mul(float128 a, float128 b, float_status_t &status)
{
int aSign, bSign, zSign;
Bit32s aExp, bExp, zExp;
Bit64u aSig0, aSig1, bSig0, bSig1, zSig0, zSig1, zSig2, zSig3;
aSig1 = extractFloat128Frac1(a);
aSig0 = extractFloat128Frac0(a);
aExp = extractFloat128Exp(a);
aSign = extractFloat128Sign(a);
bSig1 = extractFloat128Frac1(b);
bSig0 = extractFloat128Frac0(b);
bExp = extractFloat128Exp(b);
bSign = extractFloat128Sign(b);
zSign = aSign ^ bSign;
if (aExp == 0x7FFF) {
if ( (aSig0 | aSig1)
|| ((bExp == 0x7FFF) && (bSig0 | bSig1))) {
return propagateFloat128NaN(a, b, status);
}
if ((bExp | bSig0 | bSig1) == 0) goto invalid;
return packFloat128(zSign, 0x7FFF, 0, 0);
}
if (bExp == 0x7FFF) {
if (bSig0 | bSig1) return propagateFloat128NaN(a, b, status);
if ((aExp | aSig0 | aSig1) == 0) {
invalid:
float_raise(status, float_flag_invalid);
return float128_default_nan;
}
return packFloat128(zSign, 0x7FFF, 0, 0);
}
if (aExp == 0) {
if ((aSig0 | aSig1) == 0) return packFloat128(zSign, 0, 0, 0);
normalizeFloat128Subnormal(aSig0, aSig1, &aExp, &aSig0, &aSig1);
}
if (bExp == 0) {
if ((bSig0 | bSig1) == 0) return packFloat128(zSign, 0, 0, 0);
normalizeFloat128Subnormal(bSig0, bSig1, &bExp, &bSig0, &bSig1);
}
zExp = aExp + bExp - 0x4000;
aSig0 |= BX_CONST64(0x0001000000000000);
shortShift128Left(bSig0, bSig1, 16, &bSig0, &bSig1);
mul128To256(aSig0, aSig1, bSig0, bSig1, &zSig0, &zSig1, &zSig2, &zSig3);
add128(zSig0, zSig1, aSig0, aSig1, &zSig0, &zSig1);
zSig2 |= (zSig3 != 0);
if (BX_CONST64(0x0002000000000000) <= zSig0) {
shift128ExtraRightJamming(
zSig0, zSig1, zSig2, 1, &zSig0, &zSig1, &zSig2);
++zExp;
}
return roundAndPackFloat128(zSign, zExp, zSig0, zSig1, zSig2, status);
}
/*----------------------------------------------------------------------------
| Returns the result of dividing the quadruple-precision floating-point value
| `a' by the corresponding value `b'. The operation is performed according to
| the IEC/IEEE Standard for Binary Floating-Point Arithmetic.
*----------------------------------------------------------------------------*/
float128 float128_div(float128 a, float128 b, float_status_t &status)
{
int aSign, bSign, zSign;
Bit32s aExp, bExp, zExp;
Bit64u aSig0, aSig1, bSig0, bSig1, zSig0, zSig1, zSig2;
Bit64u rem0, rem1, rem2, rem3, term0, term1, term2, term3;
aSig1 = extractFloat128Frac1(a);
aSig0 = extractFloat128Frac0(a);
aExp = extractFloat128Exp(a);
aSign = extractFloat128Sign(a);
bSig1 = extractFloat128Frac1(b);
bSig0 = extractFloat128Frac0(b);
bExp = extractFloat128Exp(b);
bSign = extractFloat128Sign(b);
zSign = aSign ^ bSign;
if (aExp == 0x7FFF) {
if (aSig0 | aSig1) return propagateFloat128NaN(a, b, status);
if (bExp == 0x7FFF) {
if (bSig0 | bSig1) return propagateFloat128NaN(a, b, status);
goto invalid;
}
return packFloat128(zSign, 0x7FFF, 0, 0);
}
if (bExp == 0x7FFF) {
if (bSig0 | bSig1) return propagateFloat128NaN(a, b, status);
return packFloat128(zSign, 0, 0, 0);
}
if (bExp == 0) {
if ((bSig0 | bSig1) == 0) {
if ((aExp | aSig0 | aSig1) == 0) {
invalid:
float_raise(status, float_flag_invalid);
return float128_default_nan;
}
float_raise(status, float_flag_divbyzero);
return packFloat128(zSign, 0x7FFF, 0, 0);
}
normalizeFloat128Subnormal(bSig0, bSig1, &bExp, &bSig0, &bSig1);
}
if (aExp == 0) {
if ((aSig0 | aSig1) == 0) return packFloat128(zSign, 0, 0, 0);
normalizeFloat128Subnormal(aSig0, aSig1, &aExp, &aSig0, &aSig1);
}
zExp = aExp - bExp + 0x3FFD;
shortShift128Left(
aSig0 | BX_CONST64(0x0001000000000000), aSig1, 15, &aSig0, &aSig1);
shortShift128Left(
bSig0 | BX_CONST64(0x0001000000000000), bSig1, 15, &bSig0, &bSig1);
if (le128(bSig0, bSig1, aSig0, aSig1)) {
shift128Right(aSig0, aSig1, 1, &aSig0, &aSig1);
++zExp;
}
zSig0 = estimateDiv128To64(aSig0, aSig1, bSig0);
mul128By64To192(bSig0, bSig1, zSig0, &term0, &term1, &term2);
sub192(aSig0, aSig1, 0, term0, term1, term2, &rem0, &rem1, &rem2);
while ((Bit64s) rem0 < 0) {
--zSig0;
add192(rem0, rem1, rem2, 0, bSig0, bSig1, &rem0, &rem1, &rem2);
}
zSig1 = estimateDiv128To64(rem1, rem2, bSig0);
if ((zSig1 & 0x3FFF) <= 4) {
mul128By64To192(bSig0, bSig1, zSig1, &term1, &term2, &term3);
sub192(rem1, rem2, 0, term1, term2, term3, &rem1, &rem2, &rem3);
while ((Bit64s) rem1 < 0) {
--zSig1;
add192(rem1, rem2, rem3, 0, bSig0, bSig1, &rem1, &rem2, &rem3);
}
zSig1 |= ((rem1 | rem2 | rem3) != 0);
}
shift128ExtraRightJamming(zSig0, zSig1, 0, 15, &zSig0, &zSig1, &zSig2);
return roundAndPackFloat128(zSign, zExp, zSig0, zSig1, zSig2, status);
}
#endif
Index: Makefile.in
===================================================================
RCS file: /cvsroot/bochs/bochs/fpu/Makefile.in,v
retrieving revision 1.19.4.1
retrieving revision 1.19.4.2
diff -u -d -r1.19.4.1 -r1.19.4.2
--- Makefile.in 9 Apr 2004 12:29:49 -0000 1.19.4.1
+++ Makefile.in 22 May 2004 09:30:49 -0000 1.19.4.2
@@ -55,7 +55,8 @@
FPU_FLAGS = -DUSE_WITH_CPU_SIM $(PARANOID) $(DEBUG)
FPU_GLUE_OBJ = wmFPUemu_glue.o ferr.o fpu.o fpu_arith.o fpu_compare.o \
fpu_const.o fpu_load_store.o fpu_misc.o fpu_trans.o \
- fprem.o softfloat.o softfloatx80.o softfloat-round-pack.o
+ fprem.o softfloat.o softfloatx80.o softfloat-round-pack.o \
+ softfloat128.o poly_eval.o
# From 'C' language sources:
C_OBJS = errors.o fpu_tags.o fpu_trig.o \
@@ -163,6 +164,8 @@
../cpu/i387.h ../fpu/tag_w.h ../fpu/control_w.h fpu_proto.h
fprem.o: fprem.@...@ softfloatx80.h softfloat.h ../config.h \
softfloat-specialize.h softfloat-round-pack.h softfloat-macros.h
+poly_eval.o: poly_eval.@...@ softfloatx80.h softfloat.h ../config.h \
+ softfloat-specialize.h softfloat-round-pack.h softfloat-macros.h
ferr.o: ferr.@...@ ../bochs.h ../config.h ../osdep.h ../bx_debug/debug.h \
../bxversion.h ../gui/siminterface.h ../state_file.h ../cpu/cpu.h \
../cpu/lazy_flags.h ../cpu/hostasm.h ../cpu/apic.h ../fpu/softfloat.h \
@@ -295,6 +298,8 @@
softfloat-macros.h softfloat-specialize.h
softfloatx80.o: softfloatx80.@...@ softfloatx80.h softfloat.h ../config.h \
softfloat-specialize.h softfloat-round-pack.h softfloat-macros.h
+softfloat128.o: softfloat128.@...@ softfloatx80.h softfloat.h ../config.h \
+ softfloat-specialize.h softfloat-round-pack.h softfloat-macros.h
wmFPUemu_glue.o: wmFPUemu_glue.@...@ ../bochs.h ../config.h ../osdep.h \
../bx_debug/debug.h ../bxversion.h ../gui/siminterface.h \
../state_file.h ../cpu/cpu.h ../cpu/lazy_flags.h ../cpu/hostasm.h \
Index: div_Xsig.c
===================================================================
RCS file: /cvsroot/bochs/bochs/fpu/div_Xsig.c,v
retrieving revision 1.3.10.1
retrieving revision 1.3.10.2
diff -u -d -r1.3.10.1 -r1.3.10.2
--- div_Xsig.c 9 Apr 2004 12:29:49 -0000 1.3.10.1
+++ div_Xsig.c 22 May 2004 09:30:49 -0000 1.3.10.2
@@ -8,7 +8,6 @@
| W. Metzenthen, 22 Parker St, Ormond, Vic 3163, |
| Australia. E-mail billm@... |
| |
- | |
+---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------+
@@ -23,7 +22,6 @@
| |
| .aaaaaaaaaaaaaa / .bbbbbbbbbbbbb -> .dddddddddddd |
| |
- | |
+---------------------------------------------------------------------------*/
#include "fpu_emu.h"
@@ -178,6 +176,5 @@
result.lsw -= a.lsw;
/* Hey! we're done. */
-
*dest = result;
}
Index: softfloat-macros.h
===================================================================
RCS file: /cvsroot/bochs/bochs/fpu/Attic/softfloat-macros.h,v
retrieving revision 1.3.2.2
retrieving revision 1.3.2.3
diff -u -d -r1.3.2.2 -r1.3.2.3
--- softfloat-macros.h 11 Apr 2004 13:34:09 -0000 1.3.2.2
+++ softfloat-macros.h 22 May 2004 09:30:49 -0000 1.3.2.3
@@ -44,7 +44,7 @@
| The result is stored in the location pointed to by `zPtr'.
*----------------------------------------------------------------------------*/
-BX_CPP_INLINE void shift32RightJamming(Bit32u a, Bit16s count, Bit32u *zPtr)
+BX_CPP_INLINE void shift32RightJamming(Bit32u a, int count, Bit32u *zPtr)
{
Bit32u z;
@@ -69,20 +69,17 @@
| The result is stored in the location pointed to by `zPtr'.
*----------------------------------------------------------------------------*/
-BX_CPP_INLINE void shift64RightJamming(Bit64u a, Bit16s count, Bit64u *zPtr)
+BX_CPP_INLINE void shift64RightJamming(Bit64u a, int count, Bit64u *zPtr)
{
- Bit64u z;
-
if (count == 0) {
- z = a;
+ *zPtr = a;
}
else if (count < 64) {
- z = (a>>count) | ((a<<((-count) & 63)) != 0);
+ *zPtr = (a>>count) | ((a<<((-count) & 63)) != 0);
}
else {
- z = (a != 0);
+ *zPtr = (a != 0);
}
- *zPtr = z;
}
/*----------------------------------------------------------------------------
@@ -104,10 +101,10 @@
BX_CPP_INLINE void
shift64ExtraRightJamming(
- Bit64u a0, Bit64u a1, Bit16s count, Bit64u *z0Ptr, Bit64u *z1Ptr)
+ Bit64u a0, Bit64u a1, int count, Bit64u *z0Ptr, Bit64u *z1Ptr)
{
Bit64u z0, z1;
- Bit8s negCount = (-count) & 63;
+ int negCount = (-count) & 63;
if (count == 0) {
z1 = a1;
@@ -140,9 +137,7 @@
BX_CPP_INLINE void
add128(Bit64u a0, Bit64u a1, Bit64u b0, Bit64u b1, Bit64u *z0Ptr, Bit64u *z1Ptr)
{
- Bit64u z1;
-
- z1 = a1 + b1;
+ Bit64u z1 = a1 + b1;
*z1Ptr = z1;
*z0Ptr = a0 + b0 + (z1 < a1);
}
@@ -240,10 +235,9 @@
0x0A2D, 0x08AF, 0x075A, 0x0629, 0x051A, 0x0429, 0x0356, 0x029E,
0x0200, 0x0179, 0x0109, 0x00AF, 0x0068, 0x0034, 0x0012, 0x0002
};
- Bit8s index;
Bit32u z;
- index = (a>>27) & 15;
+ int index = (a>>27) & 15;
if (aExp & 1) {
z = 0x4000 + (a>>17) - sqrtOddAdjustments[index];
z = ((a / z)<<14) + (z<<15);
@@ -325,10 +319,10 @@
*----------------------------------------------------------------------------*/
BX_CPP_INLINE void
- shift128Right(Bit64u a0, Bit64u a1, Bit16s count, Bit64u *z0Ptr, Bit64u *z1Ptr)
+ shift128Right(Bit64u a0, Bit64u a1, int count, Bit64u *z0Ptr, Bit64u *z1Ptr)
{
Bit64u z0, z1;
- Bit8s negCount = (-count) & 63;
+ int negCount = (-count) & 63;
if (count == 0) {
z1 = a1;
@@ -359,10 +353,10 @@
BX_CPP_INLINE void
shift128RightJamming(
- Bit64u a0, Bit64u a1, Bit16s count, Bit64u *z0Ptr, Bit64u *z1Ptr)
+ Bit64u a0, Bit64u a1, int count, Bit64u *z0Ptr, Bit64u *z1Ptr)
{
Bit64u z0, z1;
- Bit8s negCount = (- count) & 63;
+ int negCount = (-count) & 63;
if (count == 0) {
z1 = a1;
@@ -397,7 +391,7 @@
BX_CPP_INLINE void
shortShift128Left(
- Bit64u a0, Bit64u a1, Bit16s count, Bit64u *z0Ptr, Bit64u *z1Ptr)
+ Bit64u a0, Bit64u a1, int count, Bit64u *z0Ptr, Bit64u *z1Ptr)
{
*z1Ptr = a1<<count;
*z0Ptr = (count == 0) ? a0 : (a0<<count) | (a1>>((-count) & 63));
@@ -510,4 +504,138 @@
#endif /* FLOATX80 */
+#ifdef FLOAT128
+
+/*----------------------------------------------------------------------------
+| Multiplies the 128-bit value formed by concatenating `a0' and `a1' by
+| `b' to obtain a 192-bit product. The product is broken into three 64-bit
+| pieces which are stored at the locations pointed to by `z0Ptr', `z1Ptr', and
+| `z2Ptr'.
+*----------------------------------------------------------------------------*/
+
+BX_CPP_INLINE void mul128By64To192(
+ Bit64u a0,
+ Bit64u a1,
+ Bit64u b,
+ Bit64u *z0Ptr,
+ Bit64u *z1Ptr,
+ Bit64u *z2Ptr
+)
+{
+ Bit64u z0, z1, z2, more1;
+
+ mul64To128(a1, b, &z1, &z2);
+ mul64To128(a0, b, &z0, &more1);
+ add128(z0, more1, 0, z1, &z0, &z1);
+ *z2Ptr = z2;
+ *z1Ptr = z1;
+ *z0Ptr = z0;
+}
+
+/*----------------------------------------------------------------------------
+| Multiplies the 128-bit value formed by concatenating `a0' and `a1' to the
+| 128-bit value formed by concatenating `b0' and `b1' to obtain a 256-bit
+| product. The product is broken into four 64-bit pieces which are stored at
+| the locations pointed to by `z0Ptr', `z1Ptr', `z2Ptr', and `z3Ptr'.
+*----------------------------------------------------------------------------*/
+
+BX_CPP_INLINE void mul128To256(
+ Bit64u a0,
+ Bit64u a1,
+ Bit64u b0,
+ Bit64u b1,
+ Bit64u *z0Ptr,
+ Bit64u *z1Ptr,
+ Bit64u *z2Ptr,
+ Bit64u *z3Ptr
+)
+{
+ Bit64u z0, z1, z2, z3;
+ Bit64u more1, more2;
+
+ mul64To128(a1, b1, &z2, &z3);
+ mul64To128(a1, b0, &z1, &more2);
+ add128(z1, more2, 0, z2, &z1, &z2);
+ mul64To128(a0, b0, &z0, &more1);
+ add128(z0, more1, 0, z1, &z0, &z1);
+ mul64To128(a0, b1, &more1, &more2);
+ add128(more1, more2, 0, z2, &more1, &z2);
+ add128(z0, z1, 0, more1, &z0, &z1);
+ *z3Ptr = z3;
+ *z2Ptr = z2;
+ *z1Ptr = z1;
+ *z0Ptr = z0;
+}
+
+
+/*----------------------------------------------------------------------------
+| Shifts the 192-bit value formed by concatenating `a0', `a1', and `a2' right
+| by 64 _plus_ the number of bits given in `count'. The shifted result is
+| at most 128 nonzero bits; these are broken into two 64-bit pieces which are
+| stored at the locations pointed to by `z0Ptr' and `z1Ptr'. The bits shifted
+| off form a third 64-bit result as follows: The _last_ bit shifted off is
+| the most-significant bit of the extra result, and the other 63 bits of the
+| extra result are all zero if and only if _all_but_the_last_ bits shifted off
+| were all zero. This extra result is stored in the location pointed to by
+| `z2Ptr'. The value of `count' can be arbitrarily large.
+| (This routine makes more sense if `a0', `a1', and `a2' are considered
+| to form a fixed-point value with binary point between `a1' and `a2'. This
+| fixed-point value is shifted right by the number of bits given in `count',
+| and the integer part of the result is returned at the locations pointed to
+| by `z0Ptr' and `z1Ptr'. The fractional part of the result may be slightly
+| corrupted as described above, and is returned at the location pointed to by
+| `z2Ptr'.)
+*----------------------------------------------------------------------------*/
+
+BX_CPP_INLINE void shift128ExtraRightJamming(
+ Bit64u a0,
+ Bit64u a1,
+ Bit64u a2,
+ int count,
+ Bit64u *z0Ptr,
+ Bit64u *z1Ptr,
+ Bit64u *z2Ptr
+)
+{
+ Bit64u z0, z1, z2;
+ int negCount = (-count) & 63;
+
+ if (count == 0) {
+ z2 = a2;
+ z1 = a1;
+ z0 = a0;
+ }
+ else {
+ if (count < 64) {
+ z2 = a1<<negCount;
+ z1 = (a0<<negCount) | (a1>>count);
+ z0 = a0>>count;
+ }
+ else {
+ if (count == 64) {
+ z2 = a1;
+ z1 = a0;
+ }
+ else {
+ a2 |= a1;
+ if (count < 128) {
+ z2 = a0<<negCount;
+ z1 = a0>>(count & 63);
+ }
+ else {
+ z2 = (count == 128) ? a0 : (a0 != 0);
+ z1 = 0;
+ }
+ }
+ z0 = 0;
+ }
+ z2 |= (a2 != 0);
+ }
+ *z2Ptr = z2;
+ *z1Ptr = z1;
+ *z0Ptr = z0;
+}
+
+#endif /* FLOAT128 */
+
#endif
Index: softfloat.h
===================================================================
RCS file: /cvsroot/bochs/bochs/fpu/Attic/softfloat.h,v
retrieving revision 1.3.2.2
retrieving revision 1.3.2.3
diff -u -d -r1.3.2.2 -r1.3.2.3
--- softfloat.h 24 Apr 2004 11:47:13 -0000 1.3.2.2
+++ softfloat.h 22 May 2004 09:30:49 -0000 1.3.2.3
@@ -295,4 +295,32 @@
#endif /* FLOATX80 */
+#ifdef FLOAT128
+
+#ifdef BX_BIG_ENDIAN
+struct float128 {
+ Bit64u hi, lo;
+};
+#else
+struct float128 {
+ Bit64u lo, hi;
+};
+#endif
+
+/*----------------------------------------------------------------------------
+| Software IEC/IEEE quadruple-precision conversion routines.
+*----------------------------------------------------------------------------*/
+float128 floatx80_to_float128(floatx80 a, float_status_t &status);
+floatx80 float128_to_floatx80(float128 a, float_status_t &status);
+
+/*----------------------------------------------------------------------------
+| Software IEC/IEEE quadruple-precision operations.
+*----------------------------------------------------------------------------*/
+float128 float128_add(float128 a, float128 b, float_status_t &status);
+float128 float128_sub(float128 a, float128 b, float_status_t &status);
+float128 float128_mul(float128 a, float128 b, float_status_t &status);
+float128 float128_div(float128 a, float128 b, float_status_t &status);
+
+#endif /* FLOAT128 */
+
#endif
|