xref: /dflybsd-src/contrib/gmp/mpf/ui_div.c (revision 86d7f5d305c6adaa56ff4582ece9859d73106103)
186d7f5d3SJohn Marino /* mpf_ui_div -- Divide an unsigned integer with a float.
286d7f5d3SJohn Marino 
386d7f5d3SJohn Marino Copyright 1993, 1994, 1995, 1996, 2000, 2001, 2002, 2004, 2005 Free Software
486d7f5d3SJohn Marino Foundation, Inc.
586d7f5d3SJohn Marino 
686d7f5d3SJohn Marino This file is part of the GNU MP Library.
786d7f5d3SJohn Marino 
886d7f5d3SJohn Marino The GNU MP Library is free software; you can redistribute it and/or modify
986d7f5d3SJohn Marino it under the terms of the GNU Lesser General Public License as published by
1086d7f5d3SJohn Marino the Free Software Foundation; either version 3 of the License, or (at your
1186d7f5d3SJohn Marino option) any later version.
1286d7f5d3SJohn Marino 
1386d7f5d3SJohn Marino The GNU MP Library is distributed in the hope that it will be useful, but
1486d7f5d3SJohn Marino WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1586d7f5d3SJohn Marino or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
1686d7f5d3SJohn Marino License for more details.
1786d7f5d3SJohn Marino 
1886d7f5d3SJohn Marino You should have received a copy of the GNU Lesser General Public License
1986d7f5d3SJohn Marino along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.  */
2086d7f5d3SJohn Marino 
2186d7f5d3SJohn Marino #include <stdio.h>  /* for NULL */
2286d7f5d3SJohn Marino #include "gmp.h"
2386d7f5d3SJohn Marino #include "gmp-impl.h"
2486d7f5d3SJohn Marino #include "longlong.h"
2586d7f5d3SJohn Marino 
2686d7f5d3SJohn Marino 
2786d7f5d3SJohn Marino void
mpf_ui_div(mpf_ptr r,unsigned long int u,mpf_srcptr v)2886d7f5d3SJohn Marino mpf_ui_div (mpf_ptr r, unsigned long int u, mpf_srcptr v)
2986d7f5d3SJohn Marino {
3086d7f5d3SJohn Marino   mp_srcptr vp;
3186d7f5d3SJohn Marino   mp_ptr rp, tp, remp, new_vp;
3286d7f5d3SJohn Marino   mp_size_t vsize;
3386d7f5d3SJohn Marino   mp_size_t rsize, prospective_rsize, zeros, tsize, high_zero;
3486d7f5d3SJohn Marino   mp_size_t sign_quotient;
3586d7f5d3SJohn Marino   mp_size_t prec;
3686d7f5d3SJohn Marino   mp_exp_t rexp;
3786d7f5d3SJohn Marino   TMP_DECL;
3886d7f5d3SJohn Marino 
3986d7f5d3SJohn Marino   vsize = v->_mp_size;
4086d7f5d3SJohn Marino   sign_quotient = vsize;
4186d7f5d3SJohn Marino   vsize = ABS (vsize);
4286d7f5d3SJohn Marino   prec = r->_mp_prec;
4386d7f5d3SJohn Marino 
4486d7f5d3SJohn Marino   if (UNLIKELY (vsize == 0))
4586d7f5d3SJohn Marino     DIVIDE_BY_ZERO;
4686d7f5d3SJohn Marino 
4786d7f5d3SJohn Marino   if (UNLIKELY (u == 0))
4886d7f5d3SJohn Marino     {
4986d7f5d3SJohn Marino       r->_mp_size = 0;
5086d7f5d3SJohn Marino       r->_mp_exp = 0;
5186d7f5d3SJohn Marino       return;
5286d7f5d3SJohn Marino     }
5386d7f5d3SJohn Marino 
5486d7f5d3SJohn Marino   TMP_MARK;
5586d7f5d3SJohn Marino   rexp = 1 - v->_mp_exp + 1;
5686d7f5d3SJohn Marino 
5786d7f5d3SJohn Marino   rp = r->_mp_d;
5886d7f5d3SJohn Marino   vp = v->_mp_d;
5986d7f5d3SJohn Marino 
6086d7f5d3SJohn Marino   prospective_rsize = 1 - vsize + 1;    /* quot from using given u,v sizes */
6186d7f5d3SJohn Marino   rsize = prec + 1;                     /* desired quot size */
6286d7f5d3SJohn Marino 
6386d7f5d3SJohn Marino   zeros = rsize - prospective_rsize;    /* padding u to give rsize */
6486d7f5d3SJohn Marino   tsize = 1 + zeros;                    /* u with zeros */
6586d7f5d3SJohn Marino 
6686d7f5d3SJohn Marino   if (WANT_TMP_DEBUG)
6786d7f5d3SJohn Marino     {
6886d7f5d3SJohn Marino       /* separate alloc blocks, for malloc debugging */
6986d7f5d3SJohn Marino       remp = TMP_ALLOC_LIMBS (vsize);
7086d7f5d3SJohn Marino       tp = TMP_ALLOC_LIMBS (tsize);
7186d7f5d3SJohn Marino       new_vp = NULL;
7286d7f5d3SJohn Marino       if (rp == vp)
7386d7f5d3SJohn Marino         new_vp = TMP_ALLOC_LIMBS (vsize);
7486d7f5d3SJohn Marino     }
7586d7f5d3SJohn Marino   else
7686d7f5d3SJohn Marino     {
7786d7f5d3SJohn Marino       /* one alloc with calculated size, for efficiency */
7886d7f5d3SJohn Marino       mp_size_t size = vsize + tsize + (rp == vp ? vsize : 0);
7986d7f5d3SJohn Marino       remp = TMP_ALLOC_LIMBS (size);
8086d7f5d3SJohn Marino       tp = remp + vsize;
8186d7f5d3SJohn Marino       new_vp = tp + tsize;
8286d7f5d3SJohn Marino     }
8386d7f5d3SJohn Marino 
8486d7f5d3SJohn Marino   /* ensure divisor doesn't overlap quotient */
8586d7f5d3SJohn Marino   if (rp == vp)
8686d7f5d3SJohn Marino     {
8786d7f5d3SJohn Marino       MPN_COPY (new_vp, vp, vsize);
8886d7f5d3SJohn Marino       vp = new_vp;
8986d7f5d3SJohn Marino     }
9086d7f5d3SJohn Marino 
9186d7f5d3SJohn Marino   MPN_ZERO (tp, tsize-1);
9286d7f5d3SJohn Marino 
9386d7f5d3SJohn Marino   tp[tsize - 1] = u & GMP_NUMB_MASK;
9486d7f5d3SJohn Marino #if BITS_PER_ULONG > GMP_NUMB_BITS
9586d7f5d3SJohn Marino   if (u > GMP_NUMB_MAX)
9686d7f5d3SJohn Marino     {
9786d7f5d3SJohn Marino       /* tsize-vsize+1 == rsize, so tsize >= rsize.  rsize == prec+1 >= 2,
9886d7f5d3SJohn Marino          so tsize >= 2, hence there's room for 2-limb u with nails */
9986d7f5d3SJohn Marino       ASSERT (tsize >= 2);
10086d7f5d3SJohn Marino       tp[tsize - 1] = u >> GMP_NUMB_BITS;
10186d7f5d3SJohn Marino       tp[tsize - 2] = u & GMP_NUMB_MASK;
10286d7f5d3SJohn Marino       rexp++;
10386d7f5d3SJohn Marino     }
10486d7f5d3SJohn Marino #endif
10586d7f5d3SJohn Marino 
10686d7f5d3SJohn Marino   ASSERT (tsize-vsize+1 == rsize);
10786d7f5d3SJohn Marino   mpn_tdiv_qr (rp, remp, (mp_size_t) 0, tp, tsize, vp, vsize);
10886d7f5d3SJohn Marino 
10986d7f5d3SJohn Marino   /* strip possible zero high limb */
11086d7f5d3SJohn Marino   high_zero = (rp[rsize-1] == 0);
11186d7f5d3SJohn Marino   rsize -= high_zero;
11286d7f5d3SJohn Marino   rexp -= high_zero;
11386d7f5d3SJohn Marino 
11486d7f5d3SJohn Marino   r->_mp_size = sign_quotient >= 0 ? rsize : -rsize;
11586d7f5d3SJohn Marino   r->_mp_exp = rexp;
11686d7f5d3SJohn Marino   TMP_FREE;
11786d7f5d3SJohn Marino }
118