fold-const.c

早期freebsd实现

   for a quotient (stored in *LQUO, *HQUO) and remainder (in *LREM, *HREM).
   CODE is a tree code for a kind of division, one of
   TRUNC_DIV_EXPR, FLOOR_DIV_EXPR, CEIL_DIV_EXPR, ROUND_DIV_EXPR
   or EXACT_DIV_EXPR
   It controls how the quotient is rounded to a integer.
   Return nonzero if the operation overflows.
   UNS nonzero says do unsigned division.  */

static int
div_and_round_double (code, uns,
		      lnum_orig, hnum_orig, lden_orig, hden_orig,
		      lquo, hquo, lrem, hrem)
     enum tree_code code;
     int uns;
     HOST_WIDE_INT lnum_orig, hnum_orig; /* num == numerator == dividend */
     HOST_WIDE_INT lden_orig, hden_orig; /* den == denominator == divisor */
     HOST_WIDE_INT *lquo, *hquo, *lrem, *hrem;
{
  int quo_neg = 0;
  short num[MAX_SHORTS + 1];	/* extra element for scaling.  */
  short den[MAX_SHORTS], quo[MAX_SHORTS];
  register int i, j, work;
  register int carry = 0;
  unsigned HOST_WIDE_INT lnum = lnum_orig;
  HOST_WIDE_INT hnum = hnum_orig;
  unsigned HOST_WIDE_INT lden = lden_orig;
  HOST_WIDE_INT hden = hden_orig;
  int overflow = 0;

  if ((hden == 0) && (lden == 0))
    abort ();

  /* calculate quotient sign and convert operands to unsigned.  */
  if (!uns) 
    {
      if (hnum < 0)
	{
	  quo_neg = ~ quo_neg;
	  /* (minimum integer) / (-1) is the only overflow case.  */
	  if (neg_double (lnum, hnum, &lnum, &hnum) && (lden & hden) == -1)
	    overflow = 1;
	}
      if (hden < 0) 
	{
	  quo_neg = ~ quo_neg;
	  neg_double (lden, hden, &lden, &hden);
	}
    }

  if (hnum == 0 && hden == 0)
    {				/* single precision */
      *hquo = *hrem = 0;
      *lquo = lnum / lden;	/* rounds toward zero since positive args */
      goto finish_up;
    }

  if (hnum == 0)
    {				/* trivial case: dividend < divisor */
      /* hden != 0 already checked.  */
      *hquo = *lquo = 0;
      *hrem = hnum;
      *lrem = lnum;
      goto finish_up;
    }

  bzero (quo, sizeof quo);

  bzero (num, sizeof num);	/* to zero 9th element */
  bzero (den, sizeof den);

  encode (num, lnum, hnum); 
  encode (den, lden, hden);

  /* This code requires more than just hden == 0.
     We also have to require that we don't need more than three bytes
     to hold CARRY.  If we ever did need four bytes to hold it, we
     would lose part of it when computing WORK on the next round.  */
  if (hden == 0 && ((lden << 8) >> 8) == lden)
    {				/* simpler algorithm */
      /* hnum != 0 already checked.  */
      for (i = MAX_SHORTS - 1; i >= 0; i--)
	{
	  work = num[i] + (carry << 8);
	  quo[i] = work / lden;
	  carry = work % lden;
	}
    }
  else {			/* full double precision,
				   with thanks to Don Knuth's
				   "Seminumerical Algorithms".  */
#define BASE 256
    int quo_est, scale, num_hi_sig, den_hi_sig, quo_hi_sig;

    /* Find the highest non-zero divisor digit.  */
    for (i = MAX_SHORTS - 1; ; i--)
      if (den[i] != 0) {
	den_hi_sig = i;
	break;
      }
    for (i = MAX_SHORTS - 1; ; i--)
      if (num[i] != 0) {
	num_hi_sig = i;
	break;
      }
    quo_hi_sig = num_hi_sig - den_hi_sig + 1;

    /* Insure that the first digit of the divisor is at least BASE/2.
       This is required by the quotient digit estimation algorithm.  */

    scale = BASE / (den[den_hi_sig] + 1);
    if (scale > 1) {		/* scale divisor and dividend */
      carry = 0;
      for (i = 0; i <= MAX_SHORTS - 1; i++) {
	work = (num[i] * scale) + carry;
	num[i] = work & 0xff;
	carry = work >> 8;
	if (num[i] != 0) num_hi_sig = i;
      }
      carry = 0;
      for (i = 0; i <= MAX_SHORTS - 1; i++) {
	work = (den[i] * scale) + carry;
	den[i] = work & 0xff;
	carry = work >> 8;
	if (den[i] != 0) den_hi_sig = i;
      }
    }

    /* Main loop */
    for (i = quo_hi_sig; i > 0; i--) {
      /* guess the next quotient digit, quo_est, by dividing the first
	 two remaining dividend digits by the high order quotient digit.
	 quo_est is never low and is at most 2 high.  */

      int num_hi;		/* index of highest remaining dividend digit */

      num_hi = i + den_hi_sig;

      work = (num[num_hi] * BASE) + (num_hi > 0 ? num[num_hi - 1] : 0);
      if (num[num_hi] != den[den_hi_sig]) {
	quo_est = work / den[den_hi_sig];
      }
      else {
	quo_est = BASE - 1;
      }

      /* refine quo_est so it's usually correct, and at most one high.   */
      while ((den[den_hi_sig - 1] * quo_est)
	     > (((work - (quo_est * den[den_hi_sig])) * BASE)
		 + ((num_hi - 1) > 0 ? num[num_hi - 2] : 0)))
	quo_est--;

      /* Try QUO_EST as the quotient digit, by multiplying the
         divisor by QUO_EST and subtracting from the remaining dividend.
	 Keep in mind that QUO_EST is the I - 1st digit.  */

      carry = 0;

      for (j = 0; j <= den_hi_sig; j++)
	{
	  int digit;

	  work = num[i + j - 1] - (quo_est * den[j]) + carry;
	  digit = work & 0xff;
	  carry = work >> 8;
	  if (digit < 0)
	    {
	      digit += BASE;
	      carry--;
	    }
	  num[i + j - 1] = digit;
	}

      /* if quo_est was high by one, then num[i] went negative and
	 we need to correct things.  */

      if (num[num_hi] < 0)
	{
	  quo_est--;
	  carry = 0;		/* add divisor back in */
	  for (j = 0; j <= den_hi_sig; j++)
	    {
	      work = num[i + j - 1] + den[j] + carry;
	      if (work > BASE)
		{
		  work -= BASE;
		  carry = 1;
		}
	      else
		{
		  carry = 0;
		}
	      num[i + j - 1] = work;
	    }
	  num [num_hi] += carry;
	}

      /* store the quotient digit.  */
      quo[i - 1] = quo_est;
    }
  }

  decode (quo, lquo, hquo);

 finish_up:
  /* if result is negative, make it so.  */
  if (quo_neg)
    neg_double (*lquo, *hquo, lquo, hquo);

  /* compute trial remainder:  rem = num - (quo * den)  */
  mul_double (*lquo, *hquo, lden_orig, hden_orig, lrem, hrem);
  neg_double (*lrem, *hrem, lrem, hrem);
  add_double (lnum_orig, hnum_orig, *lrem, *hrem, lrem, hrem);

  switch (code)
    {
    case TRUNC_DIV_EXPR:
    case TRUNC_MOD_EXPR:	/* round toward zero */
    case EXACT_DIV_EXPR:	/* for this one, it shouldn't matter */
      return overflow;

    case FLOOR_DIV_EXPR:
    case FLOOR_MOD_EXPR:	/* round toward negative infinity */
      if (quo_neg && (*lrem != 0 || *hrem != 0))   /* ratio < 0 && rem != 0 */
	{
	  /* quo = quo - 1;  */
	  add_double (*lquo, *hquo, (HOST_WIDE_INT) -1, (HOST_WIDE_INT)  -1,
		      lquo, hquo);
	}
      else return overflow;
      break;

    case CEIL_DIV_EXPR:
    case CEIL_MOD_EXPR:		/* round toward positive infinity */
      if (!quo_neg && (*lrem != 0 || *hrem != 0))  /* ratio > 0 && rem != 0 */
	{
	  add_double (*lquo, *hquo, (HOST_WIDE_INT) 1, (HOST_WIDE_INT) 0,
		      lquo, hquo);
	}
      else return overflow;
      break;
    
    case ROUND_DIV_EXPR:
    case ROUND_MOD_EXPR:	/* round to closest integer */
      {
	HOST_WIDE_INT labs_rem = *lrem, habs_rem = *hrem;
	HOST_WIDE_INT labs_den = lden, habs_den = hden, ltwice, htwice;

	/* get absolute values */
	if (*hrem < 0) neg_double (*lrem, *hrem, &labs_rem, &habs_rem);
	if (hden < 0) neg_double (lden, hden, &labs_den, &habs_den);

	/* if (2 * abs (lrem) >= abs (lden)) */
	mul_double ((HOST_WIDE_INT) 2, (HOST_WIDE_INT) 0,
		    labs_rem, habs_rem, &ltwice, &htwice);
	if (((unsigned HOST_WIDE_INT) habs_den
	     < (unsigned HOST_WIDE_INT) htwice)
	    || (((unsigned HOST_WIDE_INT) habs_den
		 == (unsigned HOST_WIDE_INT) htwice)
		&& ((HOST_WIDE_INT unsigned) labs_den
		    < (unsigned HOST_WIDE_INT) ltwice)))
	  {
	    if (*hquo < 0)
	      /* quo = quo - 1;  */
	      add_double (*lquo, *hquo,
			  (HOST_WIDE_INT) -1, (HOST_WIDE_INT) -1, lquo, hquo);
	    else
	      /* quo = quo + 1; */
	      add_double (*lquo, *hquo, (HOST_WIDE_INT) 1, (HOST_WIDE_INT) 0,
			  lquo, hquo);
	  }
	else return overflow;
      }
      break;

    default:
      abort ();
    }

  /* compute true remainder:  rem = num - (quo * den)  */
  mul_double (*lquo, *hquo, lden_orig, hden_orig, lrem, hrem);
  neg_double (*lrem, *hrem, lrem, hrem);
  add_double (lnum_orig, hnum_orig, *lrem, *hrem, lrem, hrem);
  return overflow;
}

/* Effectively truncate a real value to represent
   the nearest possible value in a narrower mode.
   The result is actually represented in the same data type as the argument,
   but its value is usually different.  */

REAL_VALUE_TYPE
real_value_truncate (mode, arg)
     enum machine_mode mode;
     REAL_VALUE_TYPE arg;
{
#ifdef __STDC__
  /* Make sure the value is actually stored in memory before we turn off
     the handler.  */
  volatile
#endif
    REAL_VALUE_TYPE value;
  jmp_buf handler, old_handler;
  int handled;

  if (setjmp (handler))
    {
      error ("floating overflow");
      return dconst0;
    }
  handled = push_float_handler (handler, old_handler);
  value = REAL_VALUE_TRUNCATE (mode, arg);
  pop_float_handler (handled, old_handler);
  return value;
}

#if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT

/* Check for infinity in an IEEE double precision number.  */

int
target_isinf (x)
     REAL_VALUE_TYPE x;
{
  /* The IEEE 64-bit double format.  */
  union {
    REAL_VALUE_TYPE d;
    struct {
      unsigned sign      :  1;
      unsigned exponent  : 11;
      unsigned mantissa1 : 20;
      unsigned mantissa2;
    } little_endian;
    struct {
      unsigned mantissa2;
      unsigned mantissa1 : 20;
      unsigned exponent  : 11;
      unsigned sign      :  1;
    } big_endian;    
  } u;

  u.d = dconstm1;
  if (u.big_endian.sign == 1)
    {
      u.d = x;
      return (u.big_endian.exponent == 2047
	      && u.big_endian.mantissa1 == 0
	      && u.big_endian.mantissa2 == 0);
    }
  else
    {
      u.d = x;
      return (u.little_endian.exponent == 2047
	      && u.little_endian.mantissa1 == 0
	      && u.little_endian.mantissa2 == 0);
    }
}

/* Check whether an IEEE double precision number is a NaN.  */

int
target_isnan (x)
     REAL_VALUE_TYPE x;
{
  /* The IEEE 64-bit double format.  */
  union {
    REAL_VALUE_TYPE d;
    struct {
      unsigned sign      :  1;
      unsigned exponent  : 11;
      unsigned mantissa1 : 20;
      unsigned mantissa2;
    } little_endian;
    struct {
      unsigned mantissa2;
      unsigned mantissa1 : 20;
      unsigned exponent  : 11;
      unsigned sign      :  1;
    } big_endian;    
  } u;

  u.d = dconstm1;
  if (u.big_endian.sign == 1)
    {
      u.d = x;
      return (u.big_endian.exponent == 2047
	      && (u.big_endian.mantissa1 != 0
		  || u.big_endian.mantissa2 != 0));
    }
  else
    {
      u.d = x;
      return (u.little_endian.exponent == 2047
	      && (u.little_endian.mantissa1 != 0
		  || u.little_endian.mantissa2 != 0));
    }
}

/* Check for a negative IEEE double precision number.  */

int
target_negative (x)
     REAL_VALUE_TYPE x;
{
  /* The IEEE 64-bit double format.  */
  union {
    REAL_VALUE_TYPE d;
    struct {
      unsigned sign      :  1;
      unsigned exponent  : 11;
      unsigned mantissa1 : 20;
      unsigned mantissa2;
    } little_endian;
    struct {
      unsigned mantissa2;
      unsigned mantissa1 : 20;
      unsigned exponent  : 11;
      unsigned sign      :  1;
    } big_endian;    
  } u;

  u.d = dconstm1;
  if (u.big_endian.sign == 1)
    {
      u.d = x;
      return u.big_endian.sign;
    }
  else
    {
      u.d = x;
      return u.little_endian.sign;
    }
}
#else /* Target not IEEE */

/* Let's assume other float formats don't have infinity.
   (This can be overridden by redefining REAL_VALUE_ISINF.)  */

target_isinf (x)
     REAL_VALUE_TYPE x;
{
  return 0;
}

/* Let's assume other float formats don't have NaNs.
   (This can be overridden by redefining REAL_VALUE_ISNAN.)  */

target_isnan (x)
     REAL_VALUE_TYPE x;
{
  return 0;
}

/* Let's assume other float formats don't have minus zero.
   (This can be overridden by redefining REAL_VALUE_NEGATIVE.)  */

target_negative (x)
     REAL_VALUE_TYPE x;
{
  return x < 0;
}
#endif /* Target not IEEE */

/* Split a tree IN into a constant and a variable part
   that could be combined with CODE to make IN.
   CODE must be a commutative arithmetic operation.
   Store the constant part into *CONP and the variable in &VARP.
   Return 1 if this was done; zero means the tree IN did not decompose
   this way.

   If CODE is PLUS_EXPR we also split trees that use MINUS_EXPR.
   Therefore, we must tell the caller whether the variable part
   was subtracted.  We do this by storing 1 or -1 into *VARSIGNP.
   The value stored is the coefficient for the variable term.
   The constant term we return should always be added;
   we negate it if necessary.  */

static int
split_tree (in, code, varp, conp, varsignp)