Annotation of OpenXM_contrib/gmp/tests/refmpf.c, Revision 1.1
1.1 ! ohara 1: /* Reference floating point routines.
! 2:
! 3: Copyright 1996, 2001 Free Software Foundation, Inc.
! 4:
! 5: This file is part of the GNU MP Library.
! 6:
! 7: The GNU MP Library is free software; you can redistribute it and/or modify
! 8: it under the terms of the GNU Lesser General Public License as published by
! 9: the Free Software Foundation; either version 2.1 of the License, or (at your
! 10: option) any later version.
! 11:
! 12: The GNU MP 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 Lesser General Public
! 15: License for more details.
! 16:
! 17: You should have received a copy of the GNU Lesser General Public License
! 18: along with the GNU MP 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 "tests.h"
! 25:
! 26:
! 27: void
! 28: refmpf_add (mpf_ptr w, mpf_srcptr u, mpf_srcptr v)
! 29: {
! 30: mp_size_t hi, lo, size;
! 31: mp_ptr ut, vt, wt;
! 32: int neg;
! 33: mp_exp_t exp;
! 34: mp_limb_t cy;
! 35: TMP_DECL (mark);
! 36:
! 37: TMP_MARK (mark);
! 38:
! 39: if (SIZ (u) == 0)
! 40: {
! 41: size = ABSIZ (v);
! 42: wt = (mp_ptr) TMP_ALLOC ((size+1) * BYTES_PER_MP_LIMB);
! 43: MPN_COPY (wt, PTR (v), size);
! 44: exp = EXP (v);
! 45: neg = SIZ (v) < 0;
! 46: goto done;
! 47: }
! 48: if (SIZ (v) == 0)
! 49: {
! 50: size = ABSIZ (u);
! 51: wt = (mp_ptr) TMP_ALLOC ((size+1) * BYTES_PER_MP_LIMB);
! 52: MPN_COPY (wt, PTR (u), size);
! 53: exp = EXP (u);
! 54: neg = SIZ (u) < 0;
! 55: goto done;
! 56: }
! 57: if ((SIZ (u) ^ SIZ (v)) < 0)
! 58: {
! 59: mpf_t tmp;
! 60: SIZ (tmp) = -SIZ (v);
! 61: EXP (tmp) = EXP (v);
! 62: PTR (tmp) = PTR (v);
! 63: refmpf_sub (w, u, tmp);
! 64: return;
! 65: }
! 66: neg = SIZ (u) < 0;
! 67:
! 68: /* Compute the significance of the hi and lo end of the result. */
! 69: hi = MAX (EXP (u), EXP (v));
! 70: lo = MIN (EXP (u) - ABSIZ (u), EXP (v) - ABSIZ (v));
! 71: size = hi - lo;
! 72: ut = (mp_ptr) TMP_ALLOC ((size + 1) * BYTES_PER_MP_LIMB);
! 73: vt = (mp_ptr) TMP_ALLOC ((size + 1) * BYTES_PER_MP_LIMB);
! 74: wt = (mp_ptr) TMP_ALLOC ((size + 1) * BYTES_PER_MP_LIMB);
! 75: MPN_ZERO (ut, size);
! 76: MPN_ZERO (vt, size);
! 77: {int off;
! 78: off = size + (EXP (u) - hi) - ABSIZ (u);
! 79: MPN_COPY (ut + off, PTR (u), ABSIZ (u));
! 80: off = size + (EXP (v) - hi) - ABSIZ (v);
! 81: MPN_COPY (vt + off, PTR (v), ABSIZ (v));
! 82: }
! 83:
! 84: cy = mpn_add_n (wt, ut, vt, size);
! 85: wt[size] = cy;
! 86: size += cy;
! 87: exp = hi + cy;
! 88:
! 89: done:
! 90: if (size > PREC (w))
! 91: {
! 92: wt += size - PREC (w);
! 93: size = PREC (w);
! 94: }
! 95: MPN_COPY (PTR (w), wt, size);
! 96: SIZ (w) = neg == 0 ? size : -size;
! 97: EXP (w) = exp;
! 98: TMP_FREE (mark);
! 99: }
! 100:
! 101: void
! 102: refmpf_sub (mpf_ptr w, mpf_srcptr u, mpf_srcptr v)
! 103: {
! 104: mp_size_t hi, lo, size;
! 105: mp_ptr ut, vt, wt;
! 106: int neg;
! 107: mp_exp_t exp;
! 108: TMP_DECL (mark);
! 109:
! 110: TMP_MARK (mark);
! 111:
! 112: if (SIZ (u) == 0)
! 113: {
! 114: size = ABSIZ (v);
! 115: wt = (mp_ptr) TMP_ALLOC ((size+1) * BYTES_PER_MP_LIMB);
! 116: MPN_COPY (wt, PTR (v), size);
! 117: exp = EXP (v);
! 118: neg = SIZ (v) > 0;
! 119: goto done;
! 120: }
! 121: if (SIZ (v) == 0)
! 122: {
! 123: size = ABSIZ (u);
! 124: wt = (mp_ptr) TMP_ALLOC ((size+1) * BYTES_PER_MP_LIMB);
! 125: MPN_COPY (wt, PTR (u), size);
! 126: exp = EXP (u);
! 127: neg = SIZ (u) < 0;
! 128: goto done;
! 129: }
! 130: if ((SIZ (u) ^ SIZ (v)) < 0)
! 131: {
! 132: mpf_t tmp;
! 133: SIZ (tmp) = -SIZ (v);
! 134: EXP (tmp) = EXP (v);
! 135: PTR (tmp) = PTR (v);
! 136: refmpf_add (w, u, tmp);
! 137: if (SIZ (u) < 0)
! 138: mpf_neg (w, w);
! 139: return;
! 140: }
! 141: neg = SIZ (u) < 0;
! 142:
! 143: /* Compute the significance of the hi and lo end of the result. */
! 144: hi = MAX (EXP (u), EXP (v));
! 145: lo = MIN (EXP (u) - ABSIZ (u), EXP (v) - ABSIZ (v));
! 146: size = hi - lo;
! 147: ut = (mp_ptr) TMP_ALLOC ((size + 1) * BYTES_PER_MP_LIMB);
! 148: vt = (mp_ptr) TMP_ALLOC ((size + 1) * BYTES_PER_MP_LIMB);
! 149: wt = (mp_ptr) TMP_ALLOC ((size + 1) * BYTES_PER_MP_LIMB);
! 150: MPN_ZERO (ut, size);
! 151: MPN_ZERO (vt, size);
! 152: {int off;
! 153: off = size + (EXP (u) - hi) - ABSIZ (u);
! 154: MPN_COPY (ut + off, PTR (u), ABSIZ (u));
! 155: off = size + (EXP (v) - hi) - ABSIZ (v);
! 156: MPN_COPY (vt + off, PTR (v), ABSIZ (v));
! 157: }
! 158:
! 159: if (mpn_cmp (ut, vt, size) >= 0)
! 160: mpn_sub_n (wt, ut, vt, size);
! 161: else
! 162: {
! 163: mpn_sub_n (wt, vt, ut, size);
! 164: neg ^= 1;
! 165: }
! 166: exp = hi;
! 167: while (size != 0 && wt[size - 1] == 0)
! 168: {
! 169: size--;
! 170: exp--;
! 171: }
! 172:
! 173: done:
! 174: if (size > PREC (w))
! 175: {
! 176: wt += size - PREC (w);
! 177: size = PREC (w);
! 178: }
! 179: MPN_COPY (PTR (w), wt, size);
! 180: SIZ (w) = neg == 0 ? size : -size;
! 181: EXP (w) = exp;
! 182: TMP_FREE (mark);
! 183: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to refmpf.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / gmp / tests