URL https://opencores.org/ocsvn/openrisc/openrisc/trunk

Subversion Repositories openrisc

[/] [openrisc/] [trunk/] [gnu-old/] [gdb-7.1/] [sim/] [rx/] [fpu.c] - Blame information for rev 842

Details | Compare with Previous | View Log


/* fpu.c --- FPU emulator for stand-alone RX simulator.
 
Copyright (C) 2008, 2009, 2010 Free Software Foundation, Inc.
Contributed by Red Hat, Inc.
 
This file is part of the GNU simulators.
 
This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 3 of the License, or
(at your option) any later version.
 
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.
 
You should have received a copy of the GNU General Public License
along with this program.  If not, see <http://www.gnu.org/licenses/>.  */
 
#include <stdio.h>
#include <stdlib.h>
 
#include "cpu.h"
#include "fpu.h"
 
/* FPU encodings are as follows:
 
   S EXPONENT MANTISSA
   1 12345678 12345678901234567890123
 
 
   1 00000000 00000000000000000000000   -0
 
   X 00000000 00000000000000000000001   Denormals
   X 00000000 11111111111111111111111
 
   X 00000001 XXXXXXXXXXXXXXXXXXXXXXX   Normals
   X 11111110 XXXXXXXXXXXXXXXXXXXXXXX
 
 
   1 11111111 00000000000000000000000   -Inf
 
   X 11111111 0XXXXXXXXXXXXXXXXXXXXXX   SNaN (X != 0)
   X 11111111 1XXXXXXXXXXXXXXXXXXXXXX   QNaN (X != 0)
 
*/
 
#define trace 0
#define tprintf if (trace) printf
 
/* Some magic numbers.  */
#define PLUS_MAX   0x7f7fffffUL
#define MINUS_MAX  0xff7fffffUL
#define PLUS_INF   0x7f800000UL
#define MINUS_INF  0xff800000UL
#define PLUS_ZERO  0x00000000UL
#define MINUS_ZERO 0x80000000UL
 
#define FP_RAISE(e) fp_raise(FPSWBITS_C##e)
static void
fp_raise (int mask)
{
  regs.r_fpsw |= mask;
  if (mask != FPSWBITS_CE)
    {
      if (regs.r_fpsw & (mask << FPSW_CESH))
        regs.r_fpsw |= (mask << FPSW_CFSH);
      if (regs.r_fpsw & FPSWBITS_FMASK)
        regs.r_fpsw |= FPSWBITS_FSUM;
      else
        regs.r_fpsw &= ~FPSWBITS_FSUM;
    }
}
 
/* We classify all numbers as one of these.  They correspond to the
   rows/colums in the exception tables.  */
typedef enum {
  FP_NORMAL,
  FP_PZERO,
  FP_NZERO,
  FP_PINFINITY,
  FP_NINFINITY,
  FP_DENORMAL,
  FP_QNAN,
  FP_SNAN
} FP_Type;
 
#if defined DEBUG0
static const char *fpt_names[] = {
  "Normal", "+0", "-0", "+Inf", "-Inf", "Denormal", "QNaN", "SNaN"
};
#endif
 
#define EXP_BIAS  127
#define EXP_ZERO -127
#define EXP_INF   128
 
#define MANT_BIAS 0x00080000UL
 
typedef struct {
  int exp;
  unsigned int mant; /* 24 bits */
  char type;
  char sign;
  fp_t orig_value;
} FP_Parts;
 
static void
fp_explode (fp_t f, FP_Parts *p)
{
  int exp, mant, sign;
 
  exp = ((f & 0x7f800000UL) >> 23);
  mant = f & 0x007fffffUL;
  sign = f & 0x80000000UL;
  /*printf("explode: %08x %x %2x %6x\n", f, sign, exp, mant);*/
 
  p->sign = sign ? -1 : 1;
  p->exp = exp - EXP_BIAS;
  p->orig_value = f;
  p->mant = mant | 0x00800000UL;
 
  if (p->exp == EXP_ZERO)
    {
      if (regs.r_fpsw & FPSWBITS_DN)
        mant = 0;
      if (mant)
        p->type = FP_DENORMAL;
      else
        {
          p->mant = 0;
          p->type = sign ? FP_NZERO : FP_PZERO;
        }
    }
  else if (p->exp == EXP_INF)
    {
      if (mant == 0)
        p->type = sign ? FP_NINFINITY : FP_PINFINITY;
      else if (mant & 0x00400000UL)
        p->type = FP_QNAN;
      else
        p->type = FP_SNAN;
    }
  else
    p->type = FP_NORMAL;
}
 
static fp_t
fp_implode (FP_Parts *p)
{
  int exp, mant;
 
  exp = p->exp + EXP_BIAS;
  mant = p->mant;
  /*printf("implode: exp %d mant 0x%x\n", exp, mant);*/
  if (p->type == FP_NORMAL)
    {
      while (mant
             && exp > 0
             && mant < 0x00800000UL)
        {
          mant <<= 1;
          exp --;
        }
      while (mant > 0x00ffffffUL)
        {
          mant >>= 1;
          exp ++;
        }
      if (exp < 0)
        {
          /* underflow */
          exp = 0;
          mant = 0;
          FP_RAISE (E);
        }
      if (exp >= 255)
        {
          /* overflow */
          exp = 255;
          mant = 0;
          FP_RAISE (O);
        }
    }
  mant &= 0x007fffffUL;
  exp &= 0xff;
  mant |= exp << 23;
  if (p->sign < 0)
    mant |= 0x80000000UL;
 
  return mant;
}
 
typedef union {
  unsigned long long ll;
  double d;
} U_d_ll;
 
static int checked_format = 0;
 
/* We assume a double format like this:
   S[1] E[11] M[52]
*/
 
static double
fp_to_double (FP_Parts *p)
{
  U_d_ll u;
 
  if (!checked_format)
    {
      u.d = 1.5;
      if (u.ll != 0x3ff8000000000000ULL)
        abort ();
      u.d = -225;
      if (u.ll != 0xc06c200000000000ULL)
        abort ();
      u.d = 10.1;
      if (u.ll != 0x4024333333333333ULL)
        abort ();
      checked_format = 1;
    }
 
  u.ll = 0;
  if (p->sign < 0)
    u.ll |= (1ULL << 63);
  /* Make sure a zero encoding stays a zero.  */
  if (p->exp != -EXP_BIAS)
    u.ll |= ((unsigned long long)p->exp + 1023ULL) << 52;
  u.ll |= (unsigned long long) (p->mant & 0x007fffffUL) << (52 - 23);
  return u.d;
}
 
static void
double_to_fp (double d, FP_Parts *p)
{
  int exp;
  U_d_ll u;
  int sign;
 
  u.d = d;
 
  sign = (u.ll & 0x8000000000000000ULL) ? 1 : 0;
  exp = u.ll >> 52;
  exp = (exp & 0x7ff);
 
  if (exp == 0)
    {
      /* A generated denormal should show up as an underflow, not
         here.  */
      if (sign)
        fp_explode (MINUS_ZERO, p);
      else
        fp_explode (PLUS_ZERO, p);
      return;
    }
 
  exp = exp - 1023;
  if ((exp + EXP_BIAS) > 254)
    {
      FP_RAISE (O);
      switch (regs.r_fpsw & FPSWBITS_RM)
        {
        case FPRM_NEAREST:
          if (sign)
            fp_explode (MINUS_INF, p);
          else
            fp_explode (PLUS_INF, p);
          break;
        case FPRM_ZERO:
          if (sign)
            fp_explode (MINUS_MAX, p);
          else
            fp_explode (PLUS_MAX, p);
          break;
        case FPRM_PINF:
          if (sign)
            fp_explode (MINUS_MAX, p);
          else
            fp_explode (PLUS_INF, p);
          break;
        case FPRM_NINF:
          if (sign)
            fp_explode (MINUS_INF, p);
          else
            fp_explode (PLUS_MAX, p);
          break;
        }
      return;
    }
  if ((exp + EXP_BIAS) < 1)
    {
      if (sign)
        fp_explode (MINUS_ZERO, p);
      else
        fp_explode (PLUS_ZERO, p);
      FP_RAISE (U);
    }
 
  p->sign = sign ? -1 : 1;
  p->exp = exp;
  p->mant = u.ll >> (52-23) & 0x007fffffUL;
  p->mant |= 0x00800000UL;
  p->type = FP_NORMAL;
 
  if (u.ll & 0x1fffffffULL)
    {
      switch (regs.r_fpsw & FPSWBITS_RM)
        {
        case FPRM_NEAREST:
          if (u.ll & 0x10000000ULL)
            p->mant ++;
          break;
        case FPRM_ZERO:
          break;
        case FPRM_PINF:
          if (sign == 1)
            p->mant ++;
          break;
        case FPRM_NINF:
          if (sign == -1)
            p->mant ++;
          break;
        }
      FP_RAISE (X);
    }
 
}
 
typedef enum {
  eNR,          /* Use the normal result.  */
  ePZ, eNZ,     /* +- zero */
  eSZ,          /* signed zero - XOR signs of ops together.  */
  eRZ,          /* +- zero depending on rounding mode.  */
  ePI, eNI,     /* +- Infinity */
  eSI,          /* signed infinity - XOR signs of ops together.  */
  eQN, eSN,     /* Quiet/Signalling NANs */
  eIn,          /* Invalid.  */
  eUn,          /* Unimplemented.  */
  eDZ,          /* Divide-by-zero.  */
  eLT,          /* less than */
  eGT,          /* greater than */
  eEQ,          /* equal to */
} FP_ExceptionCases;
 
#if defined DEBUG0
static const char *ex_names[] = {
  "NR", "PZ", "NZ", "SZ", "RZ", "PI", "NI", "SI", "QN", "SN", "IN", "Un", "DZ", "LT", "GT", "EQ"
};
#endif
 
/* This checks for all exceptional cases (not all FP exceptions) and
   returns TRUE if it is providing the result in *c.  If it returns
   FALSE, the caller should do the "normal" operation.  */
int
check_exceptions (FP_Parts *a, FP_Parts *b, fp_t *c,
                  FP_ExceptionCases ex_tab[5][5],
                  FP_ExceptionCases *case_ret)
{
  FP_ExceptionCases fpec;
 
  if (a->type == FP_SNAN
      || b->type == FP_SNAN)
    fpec = eIn;
  else if (a->type == FP_QNAN
           || b->type == FP_QNAN)
    fpec = eQN;
  else if (a->type == FP_DENORMAL
           || b->type == FP_DENORMAL)
    fpec = eUn;
  else
    fpec = ex_tab[(int)(a->type)][(int)(b->type)];
 
  /*printf("%s %s -> %s\n", fpt_names[(int)(a->type)], fpt_names[(int)(b->type)], ex_names[(int)(fpec)]);*/
 
  if (case_ret)
    *case_ret = fpec;
 
  switch (fpec)
    {
    case eNR:   /* Use the normal result.  */
      return 0;
 
    case ePZ:   /* + zero */
      *c = 0x00000000;
      return 1;
 
    case eNZ:   /* - zero */
      *c = 0x80000000;
      return 1;
 
    case eSZ:   /* signed zero */
      *c = (a->sign == b->sign) ? PLUS_ZERO : MINUS_ZERO;
      return 1;
 
    case eRZ:   /* +- zero depending on rounding mode.  */
      if ((regs.r_fpsw & FPSWBITS_RM) == FPRM_NINF)
        *c = 0x80000000;
      else
        *c = 0x00000000;
      return 1;
 
    case ePI:   /* + Infinity */
      *c = 0x7F800000;
      return 1;
 
    case eNI:   /* - Infinity */
      *c = 0xFF800000;
      return 1;
 
    case eSI:   /* sign Infinity */
      *c = (a->sign == b->sign) ? PLUS_INF : MINUS_INF;
      return 1;
 
    case eQN:   /* Quiet NANs */
      if (a->type == FP_QNAN)
        *c = a->orig_value;
      else
        *c = b->orig_value;
      return 1;
 
    case eSN:   /* Signalling NANs */
      if (a->type == FP_SNAN)
        *c = a->orig_value;
      else
        *c = b->orig_value;
      FP_RAISE (V);
      return 1;
 
    case eIn:   /* Invalid.  */
      FP_RAISE (V);
      if (a->type == FP_SNAN)
        *c = a->orig_value | 0x00400000;
      else if  (a->type == FP_SNAN)
        *c = b->orig_value | 0x00400000;
      else
        *c = 0x7fc00000;
      return 1;
 
    case eUn:   /* Unimplemented.  */
      FP_RAISE (E);
      return 1;
 
    case eDZ:   /* Division-by-zero.  */
      *c = (a->sign == b->sign) ? PLUS_INF : MINUS_INF;
      FP_RAISE (Z);
      return 1;
 
    default:
      return 0;
    }
}
 
#define CHECK_EXCEPTIONS(FPPa, FPPb, fpc, ex_tab) \
  if (check_exceptions (&FPPa, &FPPb, &fpc, ex_tab, 0))  \
    return fpc;
 
/* For each operation, we have two tables of how nonnormal cases are
   handled.  The DN=0 case is first, followed by the DN=1 case, with
   each table using the following layout: */
 
static FP_ExceptionCases ex_add_tab[5][5] = {
  /* N   +0   -0   +In  -In */
  { eNR, eNR, eNR, ePI, eNI }, /* Normal */
  { eNR, ePZ, eRZ, ePI, eNI }, /* +0   */
  { eNR, eRZ, eNZ, ePI, eNI }, /* -0   */
  { ePI, ePI, ePI, ePI, eIn }, /* +Inf */
  { eNI, eNI, eNI, eIn, eNI }, /* -Inf */
};
 
fp_t
rxfp_add (fp_t fa, fp_t fb)
{
  FP_Parts a, b, c;
  fp_t rv;
  double da, db;
 
  fp_explode (fa, &a);
  fp_explode (fb, &b);
  CHECK_EXCEPTIONS (a, b, rv, ex_add_tab);
 
  da = fp_to_double (&a);
  db = fp_to_double (&b);
  tprintf("%g + %g = %g\n", da, db, da+db);
 
  double_to_fp (da+db, &c);
  rv = fp_implode (&c);
  return rv;
}
 
static FP_ExceptionCases ex_sub_tab[5][5] = {
  /* N   +0   -0   +In  -In */
  { eNR, eNR, eNR, eNI, ePI }, /* Normal */
  { eNR, eRZ, ePZ, eNI, ePI }, /* +0   */
  { eNR, eNZ, eRZ, eNI, ePI }, /* -0   */
  { ePI, ePI, ePI, eIn, ePI }, /* +Inf */
  { eNI, eNI, eNI, eNI, eIn }, /* -Inf */
};
 
fp_t
rxfp_sub (fp_t fa, fp_t fb)
{
  FP_Parts a, b, c;
  fp_t rv;
  double da, db;
 
  fp_explode (fa, &a);
  fp_explode (fb, &b);
  CHECK_EXCEPTIONS (a, b, rv, ex_sub_tab);
 
  da = fp_to_double (&a);
  db = fp_to_double (&b);
  tprintf("%g - %g = %g\n", da, db, da-db);
 
  double_to_fp (da-db, &c);
  rv = fp_implode (&c);
 
  return rv;
}
 
static FP_ExceptionCases ex_mul_tab[5][5] = {
  /* N   +0   -0   +In  -In */
  { eNR, eNR, eNR, eSI, eSI }, /* Normal */
  { eNR, ePZ, eNZ, eIn, eIn }, /* +0   */
  { eNR, eNZ, ePZ, eIn, eIn }, /* -0   */
  { eSI, eIn, eIn, ePI, eNI }, /* +Inf */
  { eSI, eIn, eIn, eNI, ePI }, /* -Inf */
};
 
fp_t
rxfp_mul (fp_t fa, fp_t fb)
{
  FP_Parts a, b, c;
  fp_t rv;
  double da, db;
 
  fp_explode (fa, &a);
  fp_explode (fb, &b);
  CHECK_EXCEPTIONS (a, b, rv, ex_mul_tab);
 
  da = fp_to_double (&a);
  db = fp_to_double (&b);
  tprintf("%g x %g = %g\n", da, db, da*db);
 
  double_to_fp (da*db, &c);
  rv = fp_implode (&c);
 
  return rv;
}
 
static FP_ExceptionCases ex_div_tab[5][5] = {
  /* N   +0   -0   +In  -In */
  { eNR, eDZ, eDZ, eSZ, eSZ }, /* Normal */
  { eSZ, eIn, eIn, ePZ, eNZ }, /* +0   */
  { eSZ, eIn, eIn, eNZ, ePZ }, /* -0   */
  { eSI, ePI, eNI, eIn, eIn }, /* +Inf */
  { eSI, eNI, ePI, eIn, eIn }, /* -Inf */
};
 
fp_t
rxfp_div (fp_t fa, fp_t fb)
{
  FP_Parts a, b, c;
  fp_t rv;
  double da, db;
 
  fp_explode (fa, &a);
  fp_explode (fb, &b);
  CHECK_EXCEPTIONS (a, b, rv, ex_div_tab);
 
  da = fp_to_double (&a);
  db = fp_to_double (&b);
  tprintf("%g / %g = %g\n", da, db, da/db);
 
  double_to_fp (da/db, &c);
  rv = fp_implode (&c);
 
  return rv;
}
 
static FP_ExceptionCases ex_cmp_tab[5][5] = {
  /* N   +0   -0   +In  -In */
  { eNR, eNR, eNR, eLT, eGT }, /* Normal */
  { eNR, eEQ, eEQ, eLT, eGT }, /* +0   */
  { eNR, eEQ, eEQ, eLT, eGT }, /* -0   */
  { eGT, eGT, eGT, eEQ, eGT }, /* +Inf */
  { eLT, eLT, eLT, eLT, eEQ }, /* -Inf */
};
 
void
rxfp_cmp (fp_t fa, fp_t fb)
{
  FP_Parts a, b;
  fp_t c;
  FP_ExceptionCases reason;
  int flags = 0;
  double da, db;
 
  fp_explode (fa, &a);
  fp_explode (fb, &b);
 
  if (check_exceptions (&a, &b, &c, ex_cmp_tab, &reason))
    {
      if (reason == eQN)
        {
          /* Special case - incomparable.  */
          set_flags (FLAGBIT_Z | FLAGBIT_S | FLAGBIT_O, FLAGBIT_O);
          return;
        }
      return;
    }
 
  switch (reason)
    {
    case eEQ:
      flags = FLAGBIT_Z;
      break;
    case eLT:
      flags = FLAGBIT_S;
      break;
    case eGT:
      flags = 0;
      break;
    case eNR:
      da = fp_to_double (&a);
      db = fp_to_double (&b);
      tprintf("fcmp: %g cmp %g\n", da, db);
      if (da < db)
        flags = FLAGBIT_S;
      else if (da == db)
        flags = FLAGBIT_Z;
      else
        flags = 0;
      break;
    default:
      abort();
    }
 
  set_flags (FLAGBIT_Z | FLAGBIT_S | FLAGBIT_O, flags);
}
 
long
rxfp_ftoi (fp_t fa, int round_mode)
{
  FP_Parts a;
  fp_t rv;
  int sign;
  int whole_bits, frac_bits;
 
  fp_explode (fa, &a);
  sign = fa & 0x80000000UL;
 
  switch (a.type)
    {
    case FP_NORMAL:
      break;
    case FP_PZERO:
    case FP_NZERO:
      return 0;
    case FP_PINFINITY:
      FP_RAISE (V);
      return 0x7fffffffL;
    case FP_NINFINITY:
      FP_RAISE (V);
      return 0x80000000L;
    case FP_DENORMAL:
      FP_RAISE (E);
      return 0;
    case FP_QNAN:
    case FP_SNAN:
      FP_RAISE (V);
      return sign ? 0x80000000U : 0x7fffffff;
    }
 
  if (a.exp >= 31)
    {
      FP_RAISE (V);
      return sign ? 0x80000000U : 0x7fffffff;
    }
 
  a.exp -= 23;
 
  if (a.exp <= -25)
    {
      /* Less than 0.49999 */
      frac_bits = a.mant;
      whole_bits = 0;
    }
  else if (a.exp < 0)
    {
      frac_bits = a.mant << (32 + a.exp);
      whole_bits = a.mant >> (-a.exp);
    }
  else
    {
      frac_bits = 0;
      whole_bits = a.mant << a.exp;
    }
 
  if (frac_bits)
    {
      switch (round_mode & 3)
        {
        case FPRM_NEAREST:
          if (frac_bits & 0x80000000UL)
            whole_bits ++;
          break;
        case FPRM_ZERO:
          break;
        case FPRM_PINF:
          if (!sign)
            whole_bits ++;
          break;
        case FPRM_NINF:
          if (sign)
            whole_bits ++;
          break;
        }
    }
 
  rv = sign ? -whole_bits : whole_bits;
 
  return rv;
}
 
fp_t
rxfp_itof (long fa, int round_mode)
{
  fp_t rv;
  int sign = 0;
  unsigned int frac_bits;
  volatile unsigned int whole_bits;
  FP_Parts a;
 
  if (fa == 0)
    return PLUS_ZERO;
 
  if (fa < 0)
    {
      fa = -fa;
      sign = 1;
      a.sign = -1;
    }
  else
    a.sign = 1;
 
  whole_bits = fa;
  a.exp = 31;
 
  while (! (whole_bits & 0x80000000UL))
    {
      a.exp --;
      whole_bits <<= 1;
    }
  frac_bits = whole_bits & 0xff;
  whole_bits = whole_bits >> 8;
 
  if (frac_bits)
    {
      /* We must round */
      switch (round_mode & 3)
        {
        case FPRM_NEAREST:
          if (frac_bits & 0x80)
            whole_bits ++;
          break;
        case FPRM_ZERO:
          break;
        case FPRM_PINF:
          if (!sign)
            whole_bits ++;
          break;
        case FPRM_NINF:
          if (sign)
            whole_bits ++;
          break;
        }
    }
 
  a.mant = whole_bits;
  if (whole_bits & 0xff000000UL)
    {
      a.mant >>= 1;
      a.exp ++;
    }
 
  rv = fp_implode (&a);
  return rv;
}
 

Browse

Tools

Subversion Repositories openrisc

[/] [openrisc/] [trunk/] [gnu-old/] [gdb-7.1/] [sim/] [rx/] [fpu.c] - Blame information for rev 842

Line No.	Rev	Author	Line
1	227	jeremybenn	`/* fpu.c --- FPU emulator for stand-alone RX simulator.`
2
3			`Copyright (C) 2008, 2009, 2010 Free Software Foundation, Inc.`
4			`Contributed by Red Hat, Inc.`
5
6			`This file is part of the GNU simulators.`
7
8			`This program is free software; you can redistribute it and/or modify`
9			`it under the terms of the GNU General Public License as published by`
10			`the Free Software Foundation; either version 3 of the License, or`
11			`(at your option) any later version.`
12
13			`This program is distributed in the hope that it will be useful,`
14			`but WITHOUT ANY WARRANTY; without even the implied warranty of`
15			`MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the`
16			`GNU General Public License for more details.`
17
18			`You should have received a copy of the GNU General Public License`
19			`along with this program. If not, see <http://www.gnu.org/licenses/>. */`
20
21			`#include <stdio.h>`
22			`#include <stdlib.h>`
23
24			`#include "cpu.h"`
25			`#include "fpu.h"`
26
27			`/* FPU encodings are as follows:`
28
29			`S EXPONENT MANTISSA`
30			`1 12345678 12345678901234567890123`
31
32
33			`1 00000000 00000000000000000000000 -0`
34
35			`X 00000000 00000000000000000000001 Denormals`
36			`X 00000000 11111111111111111111111`
37
38			`X 00000001 XXXXXXXXXXXXXXXXXXXXXXX Normals`
39			`X 11111110 XXXXXXXXXXXXXXXXXXXXXXX`
40
41
42			`1 11111111 00000000000000000000000 -Inf`
43
44			`X 11111111 0XXXXXXXXXXXXXXXXXXXXXX SNaN (X != 0)`
45			`X 11111111 1XXXXXXXXXXXXXXXXXXXXXX QNaN (X != 0)`
46
47			`*/`
48
49			`#define trace 0`
50			`#define tprintf if (trace) printf`
51
52			`/* Some magic numbers. */`
53			`#define PLUS_MAX 0x7f7fffffUL`
54			`#define MINUS_MAX 0xff7fffffUL`
55			`#define PLUS_INF 0x7f800000UL`
56			`#define MINUS_INF 0xff800000UL`
57			`#define PLUS_ZERO 0x00000000UL`
58			`#define MINUS_ZERO 0x80000000UL`
59
60			`#define FP_RAISE(e) fp_raise(FPSWBITS_C##e)`
61			`static void`
62			`fp_raise (int mask)`
63			`{`
64			`regs.r_fpsw \|= mask;`
65			`if (mask != FPSWBITS_CE)`
66			`{`
67			`if (regs.r_fpsw & (mask << FPSW_CESH))`
68			`regs.r_fpsw \|= (mask << FPSW_CFSH);`
69			`if (regs.r_fpsw & FPSWBITS_FMASK)`
70			`regs.r_fpsw \|= FPSWBITS_FSUM;`
71			`else`
72			`regs.r_fpsw &= ~FPSWBITS_FSUM;`
73			`}`
74			`}`
75
76			`/* We classify all numbers as one of these. They correspond to the`
77			`rows/colums in the exception tables. */`
78			`typedef enum {`
79			`FP_NORMAL,`
80			`FP_PZERO,`
81			`FP_NZERO,`
82			`FP_PINFINITY,`
83			`FP_NINFINITY,`
84			`FP_DENORMAL,`
85			`FP_QNAN,`
86			`FP_SNAN`
87			`} FP_Type;`
88
89			`#if defined DEBUG0`
90			`static const char *fpt_names[] = {`
91			`"Normal", "+0", "-0", "+Inf", "-Inf", "Denormal", "QNaN", "SNaN"`
92			`};`
93			`#endif`
94
95			`#define EXP_BIAS 127`
96			`#define EXP_ZERO -127`
97			`#define EXP_INF 128`
98
99			`#define MANT_BIAS 0x00080000UL`
100
101			`typedef struct {`
102			`int exp;`
103			`unsigned int mant; /* 24 bits */`
104			`char type;`
105			`char sign;`
106			`fp_t orig_value;`
107			`} FP_Parts;`
108
109			`static void`
110			`fp_explode (fp_t f, FP_Parts *p)`
111			`{`
112			`int exp, mant, sign;`
113
114			`exp = ((f & 0x7f800000UL) >> 23);`
115			`mant = f & 0x007fffffUL;`
116			`sign = f & 0x80000000UL;`
117			`/printf("explode: %08x %x %2x %6x\n", f, sign, exp, mant);/`
118
119			`p->sign = sign ? -1 : 1;`
120			`p->exp = exp - EXP_BIAS;`
121			`p->orig_value = f;`
122			`p->mant = mant \| 0x00800000UL;`
123
124			`if (p->exp == EXP_ZERO)`
125			`{`
126			`if (regs.r_fpsw & FPSWBITS_DN)`
127			`mant = 0;`
128			`if (mant)`
129			`p->type = FP_DENORMAL;`
130			`else`
131			`{`
132			`p->mant = 0;`
133			`p->type = sign ? FP_NZERO : FP_PZERO;`
134			`}`
135			`}`
136			`else if (p->exp == EXP_INF)`
137			`{`
138			`if (mant == 0)`
139			`p->type = sign ? FP_NINFINITY : FP_PINFINITY;`
140			`else if (mant & 0x00400000UL)`
141			`p->type = FP_QNAN;`
142			`else`
143			`p->type = FP_SNAN;`
144			`}`
145			`else`
146			`p->type = FP_NORMAL;`
147			`}`
148
149			`static fp_t`
150			`fp_implode (FP_Parts *p)`
151			`{`
152			`int exp, mant;`
153
154			`exp = p->exp + EXP_BIAS;`
155			`mant = p->mant;`
156			`/printf("implode: exp %d mant 0x%x\n", exp, mant);/`
157			`if (p->type == FP_NORMAL)`
158			`{`
159			`while (mant`
160			`&& exp > 0`
161			`&& mant < 0x00800000UL)`
162			`{`
163			`mant <<= 1;`
164			`exp --;`
165			`}`
166			`while (mant > 0x00ffffffUL)`
167			`{`
168			`mant >>= 1;`
169			`exp ++;`
170			`}`
171			`if (exp < 0)`
172			`{`
173			`/* underflow */`
174			`exp = 0;`
175			`mant = 0;`
176			`FP_RAISE (E);`
177			`}`
178			`if (exp >= 255)`
179			`{`
180			`/* overflow */`
181			`exp = 255;`
182			`mant = 0;`
183			`FP_RAISE (O);`
184			`}`
185			`}`
186			`mant &= 0x007fffffUL;`
187			`exp &= 0xff;`
188			`mant \|= exp << 23;`
189			`if (p->sign < 0)`
190			`mant \|= 0x80000000UL;`
191
192			`return mant;`
193			`}`
194
195			`typedef union {`
196			`unsigned long long ll;`
197			`double d;`
198			`} U_d_ll;`
199
200			`static int checked_format = 0;`
201
202			`/* We assume a double format like this:`
203			`S[1] E[11] M[52]`
204			`*/`
205
206			`static double`
207			`fp_to_double (FP_Parts *p)`
208			`{`
209			`U_d_ll u;`
210
211			`if (!checked_format)`
212			`{`
213			`u.d = 1.5;`
214			`if (u.ll != 0x3ff8000000000000ULL)`
215			`abort ();`
216			`u.d = -225;`
217			`if (u.ll != 0xc06c200000000000ULL)`
218			`abort ();`
219			`u.d = 10.1;`
220			`if (u.ll != 0x4024333333333333ULL)`
221			`abort ();`
222			`checked_format = 1;`
223			`}`
224
225			`u.ll = 0;`
226			`if (p->sign < 0)`
227			`u.ll \|= (1ULL << 63);`
228			`/* Make sure a zero encoding stays a zero. */`
229			`if (p->exp != -EXP_BIAS)`
230			`u.ll \|= ((unsigned long long)p->exp + 1023ULL) << 52;`
231			`u.ll \|= (unsigned long long) (p->mant & 0x007fffffUL) << (52 - 23);`
232			`return u.d;`
233			`}`
234
235			`static void`
236			`double_to_fp (double d, FP_Parts *p)`
237			`{`
238			`int exp;`
239			`U_d_ll u;`
240			`int sign;`
241
242			`u.d = d;`
243
244			`sign = (u.ll & 0x8000000000000000ULL) ? 1 : 0;`
245			`exp = u.ll >> 52;`
246			`exp = (exp & 0x7ff);`
247
248			`if (exp == 0)`
249			`{`
250			`/* A generated denormal should show up as an underflow, not`
251			`here. */`
252			`if (sign)`
253			`fp_explode (MINUS_ZERO, p);`
254			`else`
255			`fp_explode (PLUS_ZERO, p);`
256			`return;`
257			`}`
258
259			`exp = exp - 1023;`
260			`if ((exp + EXP_BIAS) > 254)`
261			`{`
262			`FP_RAISE (O);`
263			`switch (regs.r_fpsw & FPSWBITS_RM)`
264			`{`
265			`case FPRM_NEAREST:`
266			`if (sign)`
267			`fp_explode (MINUS_INF, p);`
268			`else`
269			`fp_explode (PLUS_INF, p);`
270			`break;`
271			`case FPRM_ZERO:`
272			`if (sign)`
273			`fp_explode (MINUS_MAX, p);`
274			`else`
275			`fp_explode (PLUS_MAX, p);`
276			`break;`
277			`case FPRM_PINF:`
278			`if (sign)`
279			`fp_explode (MINUS_MAX, p);`
280			`else`
281			`fp_explode (PLUS_INF, p);`
282			`break;`
283			`case FPRM_NINF:`
284			`if (sign)`
285			`fp_explode (MINUS_INF, p);`
286			`else`
287			`fp_explode (PLUS_MAX, p);`
288			`break;`
289			`}`
290			`return;`
291			`}`
292			`if ((exp + EXP_BIAS) < 1)`
293			`{`
294			`if (sign)`
295			`fp_explode (MINUS_ZERO, p);`
296			`else`
297			`fp_explode (PLUS_ZERO, p);`
298			`FP_RAISE (U);`
299			`}`
300
301			`p->sign = sign ? -1 : 1;`
302			`p->exp = exp;`
303			`p->mant = u.ll >> (52-23) & 0x007fffffUL;`
304			`p->mant \|= 0x00800000UL;`
305			`p->type = FP_NORMAL;`
306
307			`if (u.ll & 0x1fffffffULL)`
308			`{`
309			`switch (regs.r_fpsw & FPSWBITS_RM)`
310			`{`
311			`case FPRM_NEAREST:`
312			`if (u.ll & 0x10000000ULL)`
313			`p->mant ++;`
314			`break;`
315			`case FPRM_ZERO:`
316			`break;`
317			`case FPRM_PINF:`
318			`if (sign == 1)`
319			`p->mant ++;`
320			`break;`
321			`case FPRM_NINF:`
322			`if (sign == -1)`
323			`p->mant ++;`
324			`break;`
325			`}`
326			`FP_RAISE (X);`
327			`}`
328
329			`}`
330
331			`typedef enum {`
332			`eNR, /* Use the normal result. */`
333			`ePZ, eNZ, /* +- zero */`
334			`eSZ, /* signed zero - XOR signs of ops together. */`
335			`eRZ, /* +- zero depending on rounding mode. */`
336			`ePI, eNI, /* +- Infinity */`
337			`eSI, /* signed infinity - XOR signs of ops together. */`
338			`eQN, eSN, /* Quiet/Signalling NANs */`
339			`eIn, /* Invalid. */`
340			`eUn, /* Unimplemented. */`
341			`eDZ, /* Divide-by-zero. */`
342			`eLT, /* less than */`
343			`eGT, /* greater than */`
344			`eEQ, /* equal to */`
345			`} FP_ExceptionCases;`
346
347			`#if defined DEBUG0`
348			`static const char *ex_names[] = {`
349			`"NR", "PZ", "NZ", "SZ", "RZ", "PI", "NI", "SI", "QN", "SN", "IN", "Un", "DZ", "LT", "GT", "EQ"`
350			`};`
351			`#endif`
352
353			`/* This checks for all exceptional cases (not all FP exceptions) and`
354			`returns TRUE if it is providing the result in *c. If it returns`
355			`FALSE, the caller should do the "normal" operation. */`
356			`int`
357			`check_exceptions (FP_Parts a, FP_Parts b, fp_t *c,`
358			`FP_ExceptionCases ex_tab[5][5],`
359			`FP_ExceptionCases *case_ret)`
360			`{`
361			`FP_ExceptionCases fpec;`
362
363			`if (a->type == FP_SNAN`
364			`\|\| b->type == FP_SNAN)`
365			`fpec = eIn;`
366			`else if (a->type == FP_QNAN`
367			`\|\| b->type == FP_QNAN)`
368			`fpec = eQN;`
369			`else if (a->type == FP_DENORMAL`
370			`\|\| b->type == FP_DENORMAL)`
371			`fpec = eUn;`
372			`else`
373			`fpec = ex_tab[(int)(a->type)][(int)(b->type)];`
374
375			`/printf("%s %s -> %s\n", fpt_names[(int)(a->type)], fpt_names[(int)(b->type)], ex_names[(int)(fpec)]);/`
376
377			`if (case_ret)`
378			`*case_ret = fpec;`
379
380			`switch (fpec)`
381			`{`
382			`case eNR: /* Use the normal result. */`
383			`return 0;`
384
385			`case ePZ: /* + zero */`
386			`*c = 0x00000000;`
387			`return 1;`
388
389			`case eNZ: /* - zero */`
390			`*c = 0x80000000;`
391			`return 1;`
392
393			`case eSZ: /* signed zero */`
394			`*c = (a->sign == b->sign) ? PLUS_ZERO : MINUS_ZERO;`
395			`return 1;`
396
397			`case eRZ: /* +- zero depending on rounding mode. */`
398			`if ((regs.r_fpsw & FPSWBITS_RM) == FPRM_NINF)`
399			`*c = 0x80000000;`
400			`else`
401			`*c = 0x00000000;`
402			`return 1;`
403
404			`case ePI: /* + Infinity */`
405			`*c = 0x7F800000;`
406			`return 1;`
407
408			`case eNI: /* - Infinity */`
409			`*c = 0xFF800000;`
410			`return 1;`
411
412			`case eSI: /* sign Infinity */`
413			`*c = (a->sign == b->sign) ? PLUS_INF : MINUS_INF;`
414			`return 1;`
415
416			`case eQN: /* Quiet NANs */`
417			`if (a->type == FP_QNAN)`
418			`*c = a->orig_value;`
419			`else`
420			`*c = b->orig_value;`
421			`return 1;`
422
423			`case eSN: /* Signalling NANs */`
424			`if (a->type == FP_SNAN)`
425			`*c = a->orig_value;`
426			`else`
427			`*c = b->orig_value;`
428			`FP_RAISE (V);`
429			`return 1;`
430
431			`case eIn: /* Invalid. */`
432			`FP_RAISE (V);`
433			`if (a->type == FP_SNAN)`
434			`*c = a->orig_value \| 0x00400000;`
435			`else if (a->type == FP_SNAN)`
436			`*c = b->orig_value \| 0x00400000;`
437			`else`
438			`*c = 0x7fc00000;`
439			`return 1;`
440
441			`case eUn: /* Unimplemented. */`
442			`FP_RAISE (E);`
443			`return 1;`
444
445			`case eDZ: /* Division-by-zero. */`
446			`*c = (a->sign == b->sign) ? PLUS_INF : MINUS_INF;`
447			`FP_RAISE (Z);`
448			`return 1;`
449
450			`default:`
451			`return 0;`
452			`}`
453			`}`
454
455			`#define CHECK_EXCEPTIONS(FPPa, FPPb, fpc, ex_tab) \`
456			`if (check_exceptions (&FPPa, &FPPb, &fpc, ex_tab, 0)) \`
457			`return fpc;`
458
459			`/* For each operation, we have two tables of how nonnormal cases are`
460			`handled. The DN=0 case is first, followed by the DN=1 case, with`
461			`each table using the following layout: */`
462
463			`static FP_ExceptionCases ex_add_tab[5][5] = {`
464			`/* N +0 -0 +In -In */`
465			`{ eNR, eNR, eNR, ePI, eNI }, /* Normal */`
466			`{ eNR, ePZ, eRZ, ePI, eNI }, /* +0 */`
467			`{ eNR, eRZ, eNZ, ePI, eNI }, /* -0 */`
468			`{ ePI, ePI, ePI, ePI, eIn }, /* +Inf */`
469			`{ eNI, eNI, eNI, eIn, eNI }, /* -Inf */`
470			`};`
471
472			`fp_t`
473			`rxfp_add (fp_t fa, fp_t fb)`
474			`{`
475			`FP_Parts a, b, c;`
476			`fp_t rv;`
477			`double da, db;`
478
479			`fp_explode (fa, &a);`
480			`fp_explode (fb, &b);`
481			`CHECK_EXCEPTIONS (a, b, rv, ex_add_tab);`
482
483			`da = fp_to_double (&a);`
484			`db = fp_to_double (&b);`
485			`tprintf("%g + %g = %g\n", da, db, da+db);`
486
487			`double_to_fp (da+db, &c);`
488			`rv = fp_implode (&c);`
489			`return rv;`
490			`}`
491
492			`static FP_ExceptionCases ex_sub_tab[5][5] = {`
493			`/* N +0 -0 +In -In */`
494			`{ eNR, eNR, eNR, eNI, ePI }, /* Normal */`
495			`{ eNR, eRZ, ePZ, eNI, ePI }, /* +0 */`
496			`{ eNR, eNZ, eRZ, eNI, ePI }, /* -0 */`
497			`{ ePI, ePI, ePI, eIn, ePI }, /* +Inf */`
498			`{ eNI, eNI, eNI, eNI, eIn }, /* -Inf */`
499			`};`
500
501			`fp_t`
502			`rxfp_sub (fp_t fa, fp_t fb)`
503			`{`
504			`FP_Parts a, b, c;`
505			`fp_t rv;`
506			`double da, db;`
507
508			`fp_explode (fa, &a);`
509			`fp_explode (fb, &b);`
510			`CHECK_EXCEPTIONS (a, b, rv, ex_sub_tab);`
511
512			`da = fp_to_double (&a);`
513			`db = fp_to_double (&b);`
514			`tprintf("%g - %g = %g\n", da, db, da-db);`
515
516			`double_to_fp (da-db, &c);`
517			`rv = fp_implode (&c);`
518
519			`return rv;`
520			`}`
521
522			`static FP_ExceptionCases ex_mul_tab[5][5] = {`
523			`/* N +0 -0 +In -In */`
524			`{ eNR, eNR, eNR, eSI, eSI }, /* Normal */`
525			`{ eNR, ePZ, eNZ, eIn, eIn }, /* +0 */`
526			`{ eNR, eNZ, ePZ, eIn, eIn }, /* -0 */`
527			`{ eSI, eIn, eIn, ePI, eNI }, /* +Inf */`
528			`{ eSI, eIn, eIn, eNI, ePI }, /* -Inf */`
529			`};`
530
531			`fp_t`
532			`rxfp_mul (fp_t fa, fp_t fb)`
533			`{`
534			`FP_Parts a, b, c;`
535			`fp_t rv;`
536			`double da, db;`
537
538			`fp_explode (fa, &a);`
539			`fp_explode (fb, &b);`
540			`CHECK_EXCEPTIONS (a, b, rv, ex_mul_tab);`
541
542			`da = fp_to_double (&a);`
543			`db = fp_to_double (&b);`
544			`tprintf("%g x %g = %g\n", da, db, da*db);`
545
546			`double_to_fp (da*db, &c);`
547			`rv = fp_implode (&c);`
548
549			`return rv;`
550			`}`
551
552			`static FP_ExceptionCases ex_div_tab[5][5] = {`
553			`/* N +0 -0 +In -In */`
554			`{ eNR, eDZ, eDZ, eSZ, eSZ }, /* Normal */`
555			`{ eSZ, eIn, eIn, ePZ, eNZ }, /* +0 */`
556			`{ eSZ, eIn, eIn, eNZ, ePZ }, /* -0 */`
557			`{ eSI, ePI, eNI, eIn, eIn }, /* +Inf */`
558			`{ eSI, eNI, ePI, eIn, eIn }, /* -Inf */`
559			`};`
560
561			`fp_t`
562			`rxfp_div (fp_t fa, fp_t fb)`
563			`{`
564			`FP_Parts a, b, c;`
565			`fp_t rv;`
566			`double da, db;`
567
568			`fp_explode (fa, &a);`
569			`fp_explode (fb, &b);`
570			`CHECK_EXCEPTIONS (a, b, rv, ex_div_tab);`
571
572			`da = fp_to_double (&a);`
573			`db = fp_to_double (&b);`
574			`tprintf("%g / %g = %g\n", da, db, da/db);`
575
576			`double_to_fp (da/db, &c);`
577			`rv = fp_implode (&c);`
578
579			`return rv;`
580			`}`
581
582			`static FP_ExceptionCases ex_cmp_tab[5][5] = {`
583			`/* N +0 -0 +In -In */`
584			`{ eNR, eNR, eNR, eLT, eGT }, /* Normal */`
585			`{ eNR, eEQ, eEQ, eLT, eGT }, /* +0 */`
586			`{ eNR, eEQ, eEQ, eLT, eGT }, /* -0 */`
587			`{ eGT, eGT, eGT, eEQ, eGT }, /* +Inf */`
588			`{ eLT, eLT, eLT, eLT, eEQ }, /* -Inf */`
589			`};`
590
591			`void`
592			`rxfp_cmp (fp_t fa, fp_t fb)`
593			`{`
594			`FP_Parts a, b;`
595			`fp_t c;`
596			`FP_ExceptionCases reason;`
597			`int flags = 0;`
598			`double da, db;`
599
600			`fp_explode (fa, &a);`
601			`fp_explode (fb, &b);`
602
603			`if (check_exceptions (&a, &b, &c, ex_cmp_tab, &reason))`
604			`{`
605			`if (reason == eQN)`
606			`{`
607			`/* Special case - incomparable. */`
608			`set_flags (FLAGBIT_Z \| FLAGBIT_S \| FLAGBIT_O, FLAGBIT_O);`
609			`return;`
610			`}`
611			`return;`
612			`}`
613
614			`switch (reason)`
615			`{`
616			`case eEQ:`
617			`flags = FLAGBIT_Z;`
618			`break;`
619			`case eLT:`
620			`flags = FLAGBIT_S;`
621			`break;`
622			`case eGT:`
623			`flags = 0;`
624			`break;`
625			`case eNR:`
626			`da = fp_to_double (&a);`
627			`db = fp_to_double (&b);`
628			`tprintf("fcmp: %g cmp %g\n", da, db);`
629			`if (da < db)`
630			`flags = FLAGBIT_S;`
631			`else if (da == db)`
632			`flags = FLAGBIT_Z;`
633			`else`
634			`flags = 0;`
635			`break;`
636			`default:`
637			`abort();`
638			`}`
639
640			`set_flags (FLAGBIT_Z \| FLAGBIT_S \| FLAGBIT_O, flags);`
641			`}`
642
643			`long`
644			`rxfp_ftoi (fp_t fa, int round_mode)`
645			`{`
646			`FP_Parts a;`
647			`fp_t rv;`
648			`int sign;`
649			`int whole_bits, frac_bits;`
650
651			`fp_explode (fa, &a);`
652			`sign = fa & 0x80000000UL;`
653
654			`switch (a.type)`
655			`{`
656			`case FP_NORMAL:`
657			`break;`
658			`case FP_PZERO:`
659			`case FP_NZERO:`
660			`return 0;`
661			`case FP_PINFINITY:`
662			`FP_RAISE (V);`
663			`return 0x7fffffffL;`
664			`case FP_NINFINITY:`
665			`FP_RAISE (V);`
666			`return 0x80000000L;`
667			`case FP_DENORMAL:`
668			`FP_RAISE (E);`
669			`return 0;`
670			`case FP_QNAN:`
671			`case FP_SNAN:`
672			`FP_RAISE (V);`
673			`return sign ? 0x80000000U : 0x7fffffff;`
674			`}`
675
676			`if (a.exp >= 31)`
677			`{`
678			`FP_RAISE (V);`
679			`return sign ? 0x80000000U : 0x7fffffff;`
680			`}`
681
682			`a.exp -= 23;`
683
684			`if (a.exp <= -25)`
685			`{`
686			`/* Less than 0.49999 */`
687			`frac_bits = a.mant;`
688			`whole_bits = 0;`
689			`}`
690			`else if (a.exp < 0)`
691			`{`
692			`frac_bits = a.mant << (32 + a.exp);`
693			`whole_bits = a.mant >> (-a.exp);`
694			`}`
695			`else`
696			`{`
697			`frac_bits = 0;`
698			`whole_bits = a.mant << a.exp;`
699			`}`
700
701			`if (frac_bits)`
702			`{`
703			`switch (round_mode & 3)`
704			`{`
705			`case FPRM_NEAREST:`
706			`if (frac_bits & 0x80000000UL)`
707			`whole_bits ++;`
708			`break;`
709			`case FPRM_ZERO:`
710			`break;`
711			`case FPRM_PINF:`
712			`if (!sign)`
713			`whole_bits ++;`
714			`break;`
715			`case FPRM_NINF:`
716			`if (sign)`
717			`whole_bits ++;`
718			`break;`
719			`}`
720			`}`
721
722			`rv = sign ? -whole_bits : whole_bits;`
723
724			`return rv;`
725			`}`
726
727			`fp_t`
728			`rxfp_itof (long fa, int round_mode)`
729			`{`
730			`fp_t rv;`
731			`int sign = 0;`
732			`unsigned int frac_bits;`
733			`volatile unsigned int whole_bits;`
734			`FP_Parts a;`
735
736			`if (fa == 0)`
737			`return PLUS_ZERO;`
738
739			`if (fa < 0)`
740			`{`
741			`fa = -fa;`
742			`sign = 1;`
743			`a.sign = -1;`
744			`}`
745			`else`
746			`a.sign = 1;`
747
748			`whole_bits = fa;`
749			`a.exp = 31;`
750
751			`while (! (whole_bits & 0x80000000UL))`
752			`{`
753			`a.exp --;`
754			`whole_bits <<= 1;`
755			`}`
756			`frac_bits = whole_bits & 0xff;`
757			`whole_bits = whole_bits >> 8;`
758
759			`if (frac_bits)`
760			`{`
761			`/* We must round */`
762			`switch (round_mode & 3)`
763			`{`
764			`case FPRM_NEAREST:`
765			`if (frac_bits & 0x80)`
766			`whole_bits ++;`
767			`break;`
768			`case FPRM_ZERO:`
769			`break;`
770			`case FPRM_PINF:`
771			`if (!sign)`
772			`whole_bits ++;`
773			`break;`
774			`case FPRM_NINF:`
775			`if (sign)`
776			`whole_bits ++;`
777			`break;`
778			`}`
779			`}`
780
781			`a.mant = whole_bits;`
782			`if (whole_bits & 0xff000000UL)`
783			`{`
784			`a.mant >>= 1;`
785			`a.exp ++;`
786			`}`
787
788			`rv = fp_implode (&a);`
789			`return rv;`
790			`}`
791