Return to set_z.c CVS log

Up to [local] / OpenXM_contrib / gmp / mpfr

Annotation of OpenXM_contrib/gmp/mpfr/set_z.c, Revision 1.1.1.1

1.1       maekawa     1: /* mpfr_set_z -- set a floating-point number from a multiple-precision integer
                      2:
                      3: Copyright (C) 1999 PolKA project, Inria Lorraine and Loria
                      4:
                      5: This file is part of the MPFR Library.
                      6:
                      7: The MPFR Library is free software; you can redistribute it and/or modify
                      8: it under the terms of the GNU Library General Public License as published by
                      9: the Free Software Foundation; either version 2 of the License, or (at your
                     10: option) any later version.
                     11:
                     12: The MPFR Library is distributed in the hope that it will be useful, but
                     13: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
                     14: or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
                     15: License for more details.
                     16:
                     17: You should have received a copy of the GNU Library General Public License
                     18: along with the MPFR Library; see the file COPYING.LIB.  If not, write to
                     19: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
                     20: MA 02111-1307, USA. */
                     21:
                     22: #include "gmp.h"
                     23: #include "gmp-impl.h"
                     24: #include "longlong.h"
                     25: #include "mpfr.h"
                     26:
                     27: /* set f to the integer z */
                     28: int
                     29: #if __STDC__
                     30: mpfr_set_z (mpfr_ptr f, mpz_srcptr z, unsigned char rnd)
                     31: #else
                     32: mpfr_set_z (f, z, rnd)
                     33:      mpfr_ptr f;
                     34:      mpz_srcptr z;
                     35:      unsigned char rnd;
                     36: #endif
                     37: {
                     38:   int fn, zn, k, dif, sign_z, sh; mp_limb_t *fp = MANT(f), *zp, cc, c2;
                     39:
                     40:   sign_z = mpz_cmp_ui(z,0);
                     41:   if (sign_z==0) return (SIZE(f)=0);
                     42:   fn = 1 + (PREC(f)-1)/BITS_PER_MP_LIMB;
                     43:   zn = ABS(SIZ(z));
                     44:   dif = zn-fn;
                     45:   zp = PTR(z);
                     46:   count_leading_zeros(k, zp[zn-1]);
                     47:   EXP(f) = zn*BITS_PER_MP_LIMB-k;
                     48:   if (SIGN(f)*sign_z<0) CHANGE_SIGN(f);
                     49:   if (dif>=0) { /* number has to be truncated */
                     50:     if (k) {
                     51:       mpn_lshift(fp, zp + dif, fn, k);
                     52:       if (dif) *fp += zp[dif-1] >> (BITS_PER_MP_LIMB-k);
                     53:     }
                     54:     else MPN_COPY(fp, zp + dif, fn);
                     55:     sh = fn*BITS_PER_MP_LIMB-PREC(f);
                     56:     cc = *fp & (((mp_limb_t)1 << sh) - 1);
                     57:     *fp = *fp & ~cc;
                     58:     if (rnd==GMP_RNDN) {
                     59:       if (sh) c2 = (mp_limb_t)1 << (sh-1);
                     60:       else { /* sh=0 */
                     61:        c2 = ((mp_limb_t)1) << (BITS_PER_MP_LIMB-1);
                     62:        dif--;
                     63:        cc = (dif>=0) ? ((zp[dif])<<k) : 0;
                     64:        if (dif>0 && k) cc += zp[dif-1] >> (BITS_PER_MP_LIMB-k);
                     65:       }
                     66:       /* now compares cc to c2 */
                     67:       if (cc>c2) { mpfr_add_one_ulp(f); return cc; }
                     68:       else if (cc<c2) goto towards_zero;
                     69:       else {
                     70:        cc=0;
                     71:        while (dif>0 && (cc=zp[dif-1])==0) dif--;
                     72:        if (cc) { mpfr_add_one_ulp(f); return cc; }
                     73:        else /* exactly in middle: inexact in both cases */
                     74:          if (*fp & ((mp_limb_t)1<<sh)) { mpfr_add_one_ulp(f); return 1; }
                     75:          else return 1;
                     76:       }
                     77:     }
                     78:     else if ((sign_z>0 && rnd==GMP_RNDU) || (sign_z<0 && rnd==GMP_RNDD)) {
                     79:       /* round towards infinity */
                     80:       /* result is exact iff all remaining bits are zero */
                     81:       if (dif>0 && cc==0) cc=zp[--dif]<<k;
                     82:       while (cc==0 && dif>0) cc=zp[--dif];
                     83:       if (cc) { mpfr_add_one_ulp(f); return 1; }
                     84:       else return 0;
                     85:     }
                     86:     else { /* round towards zero */
                     87:       /* result is exact iff all remaining bits are zero */
                     88:     towards_zero:
                     89:       if (cc==0 && dif>0) cc=zp[--dif]<<k;
                     90:       while (cc==0 && dif>0) cc=zp[--dif];
                     91:       return cc;
                     92:     }
                     93:   }
                     94:   else {
                     95:     if (k) mpn_lshift(fp-dif, zp, zn, k);
                     96:     else MPN_COPY(fp-dif, zp, zn);
                     97:     /* fill with zeroes */
                     98:     MPN_ZERO(fp, -dif);
                     99:     return 0; /* result is exact */
                    100:   }
                    101: }
                    102:
                    103:
                    104: