Blame mpf/ui_sub.c

Packit 5c3484
/* mpf_ui_sub -- Subtract a float from an unsigned long int.
Packit 5c3484
Packit 5c3484
Copyright 1993-1996, 2001, 2002, 2005, 2014 Free Software Foundation, Inc.
Packit 5c3484
Packit 5c3484
This file is part of the GNU MP Library.
Packit 5c3484
Packit 5c3484
The GNU MP Library is free software; you can redistribute it and/or modify
Packit 5c3484
it under the terms of either:
Packit 5c3484
Packit 5c3484
  * the GNU Lesser General Public License as published by the Free
Packit 5c3484
    Software Foundation; either version 3 of the License, or (at your
Packit 5c3484
    option) any later version.
Packit 5c3484
Packit 5c3484
or
Packit 5c3484
Packit 5c3484
  * the GNU General Public License as published by the Free Software
Packit 5c3484
    Foundation; either version 2 of the License, or (at your option) any
Packit 5c3484
    later version.
Packit 5c3484
Packit 5c3484
or both in parallel, as here.
Packit 5c3484
Packit 5c3484
The GNU MP Library is distributed in the hope that it will be useful, but
Packit 5c3484
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
Packit 5c3484
or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
Packit 5c3484
for more details.
Packit 5c3484
Packit 5c3484
You should have received copies of the GNU General Public License and the
Packit 5c3484
GNU Lesser General Public License along with the GNU MP Library.  If not,
Packit 5c3484
see https://www.gnu.org/licenses/.  */
Packit 5c3484
Packit 5c3484
#include "gmp.h"
Packit 5c3484
#include "gmp-impl.h"
Packit 5c3484
Packit 5c3484
void
Packit 5c3484
mpf_ui_sub (mpf_ptr r, unsigned long int u, mpf_srcptr v)
Packit 5c3484
{
Packit 5c3484
#if 1
Packit 5c3484
  __mpf_struct uu;
Packit 5c3484
  mp_limb_t ul;
Packit 5c3484
Packit 5c3484
  if (u == 0)
Packit 5c3484
    {
Packit 5c3484
      mpf_neg (r, v);
Packit 5c3484
      return;
Packit 5c3484
    }
Packit 5c3484
Packit 5c3484
  ul = u;
Packit 5c3484
  uu._mp_size = 1;
Packit 5c3484
  uu._mp_d = &ul;
Packit 5c3484
  uu._mp_exp = 1;
Packit 5c3484
  mpf_sub (r, &uu, v);
Packit 5c3484
Packit 5c3484
#else
Packit 5c3484
  mp_srcptr up, vp;
Packit 5c3484
  mp_ptr rp, tp;
Packit 5c3484
  mp_size_t usize, vsize, rsize;
Packit 5c3484
  mp_size_t prec;
Packit 5c3484
  mp_exp_t uexp;
Packit 5c3484
  mp_size_t ediff;
Packit 5c3484
  int negate;
Packit 5c3484
  mp_limb_t ulimb;
Packit 5c3484
  TMP_DECL;
Packit 5c3484
Packit 5c3484
  vsize = v->_mp_size;
Packit 5c3484
Packit 5c3484
  /* Handle special cases that don't work in generic code below.  */
Packit 5c3484
  if (u == 0)
Packit 5c3484
    {
Packit 5c3484
      mpf_neg (r, v);
Packit 5c3484
      return;
Packit 5c3484
    }
Packit 5c3484
  if (vsize == 0)
Packit 5c3484
    {
Packit 5c3484
      mpf_set_ui (r, u);
Packit 5c3484
      return;
Packit 5c3484
    }
Packit 5c3484
Packit 5c3484
  /* If signs of U and V are different, perform addition.  */
Packit 5c3484
  if (vsize < 0)
Packit 5c3484
    {
Packit 5c3484
      __mpf_struct v_negated;
Packit 5c3484
      v_negated._mp_size = -vsize;
Packit 5c3484
      v_negated._mp_exp = v->_mp_exp;
Packit 5c3484
      v_negated._mp_d = v->_mp_d;
Packit 5c3484
      mpf_add_ui (r, &v_negated, u);
Packit 5c3484
      return;
Packit 5c3484
    }
Packit 5c3484
Packit 5c3484
  /* Signs are now known to be the same.  */
Packit 5c3484
  ASSERT (vsize > 0);
Packit 5c3484
  ulimb = u;
Packit 5c3484
  /* Make U be the operand with the largest exponent.  */
Packit 5c3484
  negate = 1 < v->_mp_exp;
Packit 5c3484
  prec = r->_mp_prec + negate;
Packit 5c3484
  rp = r->_mp_d;
Packit 5c3484
  if (negate)
Packit 5c3484
    {
Packit 5c3484
      usize = vsize;
Packit 5c3484
      vsize = 1;
Packit 5c3484
      up = v->_mp_d;
Packit 5c3484
      vp = &ulimb;
Packit 5c3484
      uexp = v->_mp_exp;
Packit 5c3484
      ediff = uexp - 1;
Packit 5c3484
Packit 5c3484
      /* If U extends beyond PREC, ignore the part that does.  */
Packit 5c3484
      if (usize > prec)
Packit 5c3484
	{
Packit 5c3484
	  up += usize - prec;
Packit 5c3484
	  usize = prec;
Packit 5c3484
	}
Packit 5c3484
      ASSERT (ediff > 0);
Packit 5c3484
    }
Packit 5c3484
  else
Packit 5c3484
    {
Packit 5c3484
      vp = v->_mp_d;
Packit 5c3484
      ediff = 1 - v->_mp_exp;
Packit 5c3484
  /* Ignore leading limbs in U and V that are equal.  Doing
Packit 5c3484
     this helps increase the precision of the result.  */
Packit 5c3484
      if (ediff == 0 && ulimb == vp[vsize - 1])
Packit 5c3484
	{
Packit 5c3484
	  usize = 0;
Packit 5c3484
	  vsize--;
Packit 5c3484
	  uexp = 0;
Packit 5c3484
	  /* Note that V might now have leading zero limbs.
Packit 5c3484
	     In that case we have to adjust uexp.  */
Packit 5c3484
	  for (;;)
Packit 5c3484
	    {
Packit 5c3484
	      if (vsize == 0) {
Packit 5c3484
		rsize = 0;
Packit 5c3484
		uexp = 0;
Packit 5c3484
		goto done;
Packit 5c3484
	      }
Packit 5c3484
	      if ( vp[vsize - 1] != 0)
Packit 5c3484
		break;
Packit 5c3484
	      vsize--, uexp--;
Packit 5c3484
	    }
Packit 5c3484
	}
Packit 5c3484
      else
Packit 5c3484
	{
Packit 5c3484
	  usize = 1;
Packit 5c3484
	  uexp = 1;
Packit 5c3484
	  up = &ulimb;
Packit 5c3484
	}
Packit 5c3484
      ASSERT (usize <= prec);
Packit 5c3484
    }
Packit 5c3484
Packit 5c3484
  if (ediff >= prec)
Packit 5c3484
    {
Packit 5c3484
      /* V completely cancelled.  */
Packit 5c3484
      if (rp != up)
Packit 5c3484
	MPN_COPY (rp, up, usize);
Packit 5c3484
      rsize = usize;
Packit 5c3484
    }
Packit 5c3484
  else
Packit 5c3484
    {
Packit 5c3484
  /* If V extends beyond PREC, ignore the part that does.
Packit 5c3484
     Note that this can make vsize neither zero nor negative.  */
Packit 5c3484
  if (vsize + ediff > prec)
Packit 5c3484
    {
Packit 5c3484
      vp += vsize + ediff - prec;
Packit 5c3484
      vsize = prec - ediff;
Packit 5c3484
    }
Packit 5c3484
Packit 5c3484
      /* Locate the least significant non-zero limb in (the needed
Packit 5c3484
	 parts of) U and V, to simplify the code below.  */
Packit 5c3484
      ASSERT (vsize > 0);
Packit 5c3484
      for (;;)
Packit 5c3484
	{
Packit 5c3484
	  if (vp[0] != 0)
Packit 5c3484
	    break;
Packit 5c3484
	  vp++, vsize--;
Packit 5c3484
	  if (vsize == 0)
Packit 5c3484
	    {
Packit 5c3484
	      MPN_COPY (rp, up, usize);
Packit 5c3484
	      rsize = usize;
Packit 5c3484
	      goto done;
Packit 5c3484
	    }
Packit 5c3484
	}
Packit 5c3484
      for (;;)
Packit 5c3484
	{
Packit 5c3484
	  if (usize == 0)
Packit 5c3484
	    {
Packit 5c3484
	      MPN_COPY (rp, vp, vsize);
Packit 5c3484
	      rsize = vsize;
Packit 5c3484
	      negate ^= 1;
Packit 5c3484
	      goto done;
Packit 5c3484
	    }
Packit 5c3484
	  if (up[0] != 0)
Packit 5c3484
	    break;
Packit 5c3484
	  up++, usize--;
Packit 5c3484
	}
Packit 5c3484
Packit 5c3484
      ASSERT (usize > 0 && vsize > 0);
Packit 5c3484
      TMP_MARK;
Packit 5c3484
Packit 5c3484
      tp = TMP_ALLOC_LIMBS (prec);
Packit 5c3484
Packit 5c3484
      /* uuuu     |  uuuu     |  uuuu     |  uuuu     |  uuuu    */
Packit 5c3484
      /* vvvvvvv  |  vv       |    vvvvv  |    v      |       vv */
Packit 5c3484
Packit 5c3484
      if (usize > ediff)
Packit 5c3484
	{
Packit 5c3484
	  /* U and V partially overlaps.  */
Packit 5c3484
	  if (ediff == 0)
Packit 5c3484
	    {
Packit 5c3484
	      ASSERT (usize == 1 && vsize >= 1 && ulimb == *up); /* usize is 1>ediff, vsize >= 1 */
Packit 5c3484
	      if (1 < vsize)
Packit 5c3484
		{
Packit 5c3484
		  /* u        */
Packit 5c3484
		  /* vvvvvvv  */
Packit 5c3484
		  rsize = vsize;
Packit 5c3484
		  vsize -= 1;
Packit 5c3484
		  /* mpn_cmp (up, vp + vsize - usize, usize) > 0 */
Packit 5c3484
		  if (ulimb > vp[vsize])
Packit 5c3484
		    {
Packit 5c3484
		      tp[vsize] = ulimb - vp[vsize] - 1;
Packit 5c3484
		      ASSERT_CARRY (mpn_neg (tp, vp, vsize));
Packit 5c3484
		    }
Packit 5c3484
		  else
Packit 5c3484
		    {
Packit 5c3484
		      /* vvvvvvv  */  /* Swap U and V. */
Packit 5c3484
		      /* u        */
Packit 5c3484
		      MPN_COPY (tp, vp, vsize);
Packit 5c3484
		      tp[vsize] = vp[vsize] - ulimb;
Packit 5c3484
		      negate = 1;
Packit 5c3484
		    }
Packit 5c3484
		}
Packit 5c3484
	      else /* vsize == usize == 1 */
Packit 5c3484
		{
Packit 5c3484
		  /* u     */
Packit 5c3484
		  /* v     */
Packit 5c3484
		  rsize = 1;
Packit 5c3484
		  negate = ulimb < vp[0];
Packit 5c3484
		  tp[0] = negate ? vp[0] - ulimb: ulimb - vp[0];
Packit 5c3484
		}
Packit 5c3484
	    }
Packit 5c3484
	  else
Packit 5c3484
	    {
Packit 5c3484
	      ASSERT (vsize + ediff <= usize);
Packit 5c3484
	      ASSERT (vsize == 1 && usize >= 2 && ulimb == *vp);
Packit 5c3484
		{
Packit 5c3484
		  /* uuuu     */
Packit 5c3484
		  /*   v      */
Packit 5c3484
		  mp_size_t size;
Packit 5c3484
		  size = usize - ediff - 1;
Packit 5c3484
		  MPN_COPY (tp, up, size);
Packit 5c3484
		  ASSERT_NOCARRY (mpn_sub_1 (tp + size, up + size, usize - size, ulimb));
Packit 5c3484
		  rsize = usize;
Packit 5c3484
		}
Packit 5c3484
		/* Other cases are not possible */
Packit 5c3484
		/* uuuu     */
Packit 5c3484
		/*   vvvvv  */
Packit 5c3484
	    }
Packit 5c3484
	}
Packit 5c3484
      else
Packit 5c3484
	{
Packit 5c3484
	  /* uuuu     */
Packit 5c3484
	  /*      vv  */
Packit 5c3484
	  mp_size_t size, i;
Packit 5c3484
	  ASSERT_CARRY (mpn_neg (tp, vp, vsize));
Packit 5c3484
	  rsize = vsize + ediff;
Packit 5c3484
	  size = rsize - usize;
Packit 5c3484
	  for (i = vsize; i < size; i++)
Packit 5c3484
	    tp[i] = GMP_NUMB_MAX;
Packit 5c3484
	  ASSERT_NOCARRY (mpn_sub_1 (tp + size, up, usize, CNST_LIMB (1)));
Packit 5c3484
	}
Packit 5c3484
Packit 5c3484
      /* Full normalize.  Optimize later.  */
Packit 5c3484
      while (rsize != 0 && tp[rsize - 1] == 0)
Packit 5c3484
	{
Packit 5c3484
	  rsize--;
Packit 5c3484
	  uexp--;
Packit 5c3484
	}
Packit 5c3484
      MPN_COPY (rp, tp, rsize);
Packit 5c3484
      TMP_FREE;
Packit 5c3484
    }
Packit 5c3484
Packit 5c3484
 done:
Packit 5c3484
  r->_mp_size = negate ? -rsize : rsize;
Packit 5c3484
  r->_mp_exp = uexp;
Packit 5c3484
#endif
Packit 5c3484
}