yann@1250
|
1 |
Original patch from: mpf_eq.diff
|
yann@1250
|
2 |
|
yann@1250
|
3 |
-= BEGIN original header =-
|
yann@1250
|
4 |
-= END original header =-
|
yann@1250
|
5 |
|
yann@1250
|
6 |
diff -durN gmp-4.2.4.orig/doc/gmp.texi gmp-4.2.4/doc/gmp.texi
|
yann@1250
|
7 |
--- gmp-4.2.4.orig/doc/gmp.texi 2008-09-18 17:36:14.000000000 +0200
|
yann@1250
|
8 |
+++ gmp-4.2.4/doc/gmp.texi 2009-03-08 18:36:16.000000000 +0100
|
yann@1250
|
9 |
@@ -4849,9 +4849,12 @@
|
yann@1250
|
10 |
equal, zero otherwise. I.e., test if @var{op1} and @var{op2} are approximately
|
yann@1250
|
11 |
equal.
|
yann@1250
|
12 |
|
yann@1250
|
13 |
-Caution: Currently only whole limbs are compared, and only in an exact
|
yann@1250
|
14 |
-fashion. In the future values like 1000 and 0111 may be considered the same
|
yann@1250
|
15 |
-to 3 bits (on the basis that their difference is that small).
|
yann@1250
|
16 |
+Caution 1: All version of GMP up to version 4.2.4 compared just whole limbs,
|
yann@1250
|
17 |
+meaning sometimes more than @var{op3} bits, sometimes fewer.
|
yann@1250
|
18 |
+
|
yann@1250
|
19 |
+Caution 2: This function will consider XXX11...111 and XX100...000 different,
|
yann@1250
|
20 |
+even if ... is replaced by a semi-infinite number of bits. Such numbers are
|
yann@1250
|
21 |
+really just one ulp off, and should be considered equal.
|
yann@1250
|
22 |
@end deftypefun
|
yann@1250
|
23 |
|
yann@1250
|
24 |
@deftypefun void mpf_reldiff (mpf_t @var{rop}, mpf_t @var{op1}, mpf_t @var{op2})
|
yann@1250
|
25 |
diff -durN gmp-4.2.4.orig/mpf/eq.c gmp-4.2.4/mpf/eq.c
|
yann@1250
|
26 |
--- gmp-4.2.4.orig/mpf/eq.c 2007-08-30 20:31:40.000000000 +0200
|
yann@1250
|
27 |
+++ gmp-4.2.4/mpf/eq.c 2009-03-08 18:36:16.000000000 +0100
|
yann@1250
|
28 |
@@ -1,6 +1,6 @@
|
yann@1250
|
29 |
/* mpf_eq -- Compare two floats up to a specified bit #.
|
yann@1250
|
30 |
|
yann@1250
|
31 |
-Copyright 1993, 1995, 1996, 2001, 2002 Free Software Foundation, Inc.
|
yann@1250
|
32 |
+Copyright 1993, 1995, 1996, 2001, 2002, 2008 Free Software Foundation, Inc.
|
yann@1250
|
33 |
|
yann@1250
|
34 |
This file is part of the GNU MP Library.
|
yann@1250
|
35 |
|
yann@1250
|
36 |
@@ -19,6 +19,7 @@
|
yann@1250
|
37 |
|
yann@1250
|
38 |
#include "gmp.h"
|
yann@1250
|
39 |
#include "gmp-impl.h"
|
yann@1250
|
40 |
+#include "longlong.h"
|
yann@1250
|
41 |
|
yann@1250
|
42 |
int
|
yann@1250
|
43 |
mpf_eq (mpf_srcptr u, mpf_srcptr v, unsigned long int n_bits)
|
yann@1250
|
44 |
@@ -26,6 +27,8 @@
|
yann@1250
|
45 |
mp_srcptr up, vp;
|
yann@1250
|
46 |
mp_size_t usize, vsize, size, i;
|
yann@1250
|
47 |
mp_exp_t uexp, vexp;
|
yann@1250
|
48 |
+ mp_limb_t diff;
|
yann@1250
|
49 |
+ int cnt;
|
yann@1250
|
50 |
|
yann@1250
|
51 |
uexp = u->_mp_exp;
|
yann@1250
|
52 |
vexp = v->_mp_exp;
|
yann@1250
|
53 |
@@ -53,10 +56,8 @@
|
yann@1250
|
54 |
/* U and V have the same sign and are both non-zero. */
|
yann@1250
|
55 |
|
yann@1250
|
56 |
/* 2. Are the exponents different? */
|
yann@1250
|
57 |
- if (uexp > vexp)
|
yann@1250
|
58 |
- return 0; /* ??? handle (uexp = vexp + 1) */
|
yann@1250
|
59 |
- if (vexp > uexp)
|
yann@1250
|
60 |
- return 0; /* ??? handle (vexp = uexp + 1) */
|
yann@1250
|
61 |
+ if (uexp != vexp)
|
yann@1250
|
62 |
+ return 0;
|
yann@1250
|
63 |
|
yann@1250
|
64 |
usize = ABS (usize);
|
yann@1250
|
65 |
vsize = ABS (vsize);
|
yann@1250
|
66 |
@@ -93,17 +94,26 @@
|
yann@1250
|
67 |
size = usize;
|
yann@1250
|
68 |
}
|
yann@1250
|
69 |
|
yann@1250
|
70 |
- if (size > (n_bits + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS)
|
yann@1250
|
71 |
- size = (n_bits + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS;
|
yann@1250
|
72 |
+ up += usize; /* point just above most significant limb */
|
yann@1250
|
73 |
+ vp += vsize; /* point just above most significant limb */
|
yann@1250
|
74 |
|
yann@1250
|
75 |
- up += usize - size;
|
yann@1250
|
76 |
- vp += vsize - size;
|
yann@1250
|
77 |
+ count_leading_zeros (cnt, up[-1]);
|
yann@1250
|
78 |
+ if ((vp[-1] >> (GMP_LIMB_BITS - 1 - cnt)) != 1)
|
yann@1250
|
79 |
+ return 0; /* msb positions different */
|
yann@1250
|
80 |
|
yann@1250
|
81 |
- for (i = size - 1; i >= 0; i--)
|
yann@1250
|
82 |
+ n_bits += cnt - GMP_NAIL_BITS;
|
yann@1250
|
83 |
+
|
yann@1250
|
84 |
+ size = MIN (size, (n_bits + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS);
|
yann@1250
|
85 |
+
|
yann@1250
|
86 |
+ up -= size; /* point at least significant relevant limb */
|
yann@1250
|
87 |
+ vp -= size; /* point at least significant relevant limb */
|
yann@1250
|
88 |
+
|
yann@1250
|
89 |
+ for (i = size - 1; i > 0; i--)
|
yann@1250
|
90 |
{
|
yann@1250
|
91 |
if (up[i] != vp[i])
|
yann@1250
|
92 |
return 0;
|
yann@1250
|
93 |
}
|
yann@1250
|
94 |
|
yann@1250
|
95 |
- return 1;
|
yann@1250
|
96 |
+ diff = (up[0] ^ vp[0]) >> GMP_NUMB_BITS - 1 - (n_bits - 1) % GMP_NUMB_BITS;
|
yann@1250
|
97 |
+ return diff == 0;
|
yann@1250
|
98 |
}
|
yann@1250
|
99 |
diff -durN gmp-4.2.4.orig/tests/cxx/t-prec.cc gmp-4.2.4/tests/cxx/t-prec.cc
|
yann@1250
|
100 |
--- gmp-4.2.4.orig/tests/cxx/t-prec.cc 2007-09-01 12:09:03.000000000 +0200
|
yann@1250
|
101 |
+++ gmp-4.2.4/tests/cxx/t-prec.cc 2009-03-08 18:36:16.000000000 +0100
|
yann@1250
|
102 |
@@ -1,6 +1,6 @@
|
yann@1250
|
103 |
/* Test precision of mpf_class expressions.
|
yann@1250
|
104 |
|
yann@1250
|
105 |
-Copyright 2001, 2002, 2003 Free Software Foundation, Inc.
|
yann@1250
|
106 |
+Copyright 2001, 2002, 2003, 2008 Free Software Foundation, Inc.
|
yann@1250
|
107 |
|
yann@1250
|
108 |
This file is part of the GNU MP Library.
|
yann@1250
|
109 |
|
yann@1250
|
110 |
@@ -61,7 +61,7 @@
|
yann@1250
|
111 |
g = 1 / f;
|
yann@1250
|
112 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
113 |
(g, "0.11111 11111 11111 11111 11111 11111 11111 11111 11111 11111"
|
yann@1250
|
114 |
- " 11111 11111 11111 11111 11111 11", very_large_prec);
|
yann@1250
|
115 |
+ " 11111 11111 11111 11111 11111 111", very_large_prec);
|
yann@1250
|
116 |
}
|
yann@1250
|
117 |
{
|
yann@1250
|
118 |
mpf_class f(15.0, large_prec);
|
yann@1250
|
119 |
@@ -69,7 +69,7 @@
|
yann@1250
|
120 |
g = 1 / f;
|
yann@1250
|
121 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
122 |
(g, "0.06666 66666 66666 66666 66666 66666 66666 66666 66666 66666"
|
yann@1250
|
123 |
- " 66666 66666 66666 66666 66666 67", very_large_prec);
|
yann@1250
|
124 |
+ " 66666 66666 66666 66666 66666 667", very_large_prec);
|
yann@1250
|
125 |
}
|
yann@1250
|
126 |
|
yann@1250
|
127 |
// compound expressions
|
yann@1250
|
128 |
@@ -94,14 +94,14 @@
|
yann@1250
|
129 |
i = f / g + h;
|
yann@1250
|
130 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
131 |
(i, "15.33333 33333 33333 33333 33333 33333 33333 33333 33333 33333"
|
yann@1250
|
132 |
- " 33333 33333 33333 333", very_large_prec);
|
yann@1250
|
133 |
+ " 33333 33333 33333 33333 33333 3", very_large_prec);
|
yann@1250
|
134 |
}
|
yann@1250
|
135 |
{
|
yann@1250
|
136 |
mpf_class f(3.0, small_prec);
|
yann@1250
|
137 |
mpf_class g(-(1 + f) / 3, very_large_prec);
|
yann@1250
|
138 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
139 |
(g, "-1.33333 33333 33333 33333 33333 33333 33333 33333 33333 33333"
|
yann@1250
|
140 |
- " 33333 33333 33333 333", very_large_prec);
|
yann@1250
|
141 |
+ " 33333 33333 33333 33333 33333 33", very_large_prec);
|
yann@1250
|
142 |
}
|
yann@1250
|
143 |
{
|
yann@1250
|
144 |
mpf_class f(9.0, medium_prec);
|
yann@1250
|
145 |
@@ -117,7 +117,7 @@
|
yann@1250
|
146 |
g = hypot(1 + 5 / f, 1.0);
|
yann@1250
|
147 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
148 |
(g, "1.66666 66666 66666 66666 66666 66666 66666 66666 66666 66666"
|
yann@1250
|
149 |
- " 66666 66666 66666 667", very_large_prec);
|
yann@1250
|
150 |
+ " 66666 66666 66666 66666 66666 67", very_large_prec);
|
yann@1250
|
151 |
}
|
yann@1250
|
152 |
|
yann@1250
|
153 |
// compound assignments
|
yann@1250
|
154 |
@@ -142,7 +142,7 @@
|
yann@1250
|
155 |
mpf_class g(0.0, very_large_prec);
|
yann@1250
|
156 |
g = mpf_class(1 / f);
|
yann@1250
|
157 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
158 |
- (g, "0.11111 11111 11111 11111 11111 11111 11111 111", medium_prec);
|
yann@1250
|
159 |
+ (g, "0.11111 11111 11111 11111 11111 11111 11111 1111", medium_prec);
|
yann@1250
|
160 |
}
|
yann@1250
|
161 |
{
|
yann@1250
|
162 |
mpf_class f(15.0, large_prec);
|
yann@1250
|
163 |
@@ -150,7 +150,7 @@
|
yann@1250
|
164 |
g = mpf_class(1 / f);
|
yann@1250
|
165 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
166 |
(g, "0.06666 66666 66666 66666 66666 66666 66666 66666 66666 66666"
|
yann@1250
|
167 |
- " 66666 667", large_prec);
|
yann@1250
|
168 |
+ " 66666 6667", large_prec);
|
yann@1250
|
169 |
}
|
yann@1250
|
170 |
|
yann@1250
|
171 |
{
|
yann@1250
|
172 |
@@ -158,7 +158,8 @@
|
yann@1250
|
173 |
mpf_class h(0.0, very_large_prec);
|
yann@1250
|
174 |
h = mpf_class(f / g + 1, large_prec);
|
yann@1250
|
175 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
176 |
- (h, "1.33333 33333 33333 33333 33333 33333 33333 33333 33333 3333",
|
yann@1250
|
177 |
+ (h, "1.33333 33333 33333 33333 33333 33333 33333 33333 33333 33333"
|
yann@1250
|
178 |
+ " 33333 333",
|
yann@1250
|
179 |
large_prec);
|
yann@1250
|
180 |
}
|
yann@1250
|
181 |
|
yann@1250
|
182 |
@@ -170,7 +171,7 @@
|
yann@1250
|
183 |
g = f - q;
|
yann@1250
|
184 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
185 |
(g, "2.66666 66666 66666 66666 66666 66666 66666 66666 66666 66666"
|
yann@1250
|
186 |
- " 66666 66666 66666 667", very_large_prec);
|
yann@1250
|
187 |
+ " 66666 66666 66666 66666 66666 67", very_large_prec);
|
yann@1250
|
188 |
}
|
yann@1250
|
189 |
|
yann@1250
|
190 |
{
|
yann@1250
|
191 |
@@ -179,7 +180,8 @@
|
yann@1250
|
192 |
mpf_class g(0.0, very_large_prec);
|
yann@1250
|
193 |
g = mpf_class(f - q, large_prec);
|
yann@1250
|
194 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
195 |
- (g, "2.66666 66666 66666 66666 66666 66666 66666 66666 66666 6667",
|
yann@1250
|
196 |
+ (g, "2.66666 66666 66666 66666 66666 66666 66666 66666 66666 66666"
|
yann@1250
|
197 |
+ " 66666 667",
|
yann@1250
|
198 |
large_prec);
|
yann@1250
|
199 |
}
|
yann@1250
|
200 |
{
|
yann@1250
|
201 |
@@ -188,7 +190,7 @@
|
yann@1250
|
202 |
mpf_class g(0.0, very_large_prec);
|
yann@1250
|
203 |
g = mpf_class(f - q);
|
yann@1250
|
204 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
205 |
- (g, "2.66666 66666 66666 66666 66666 6667", medium_prec);
|
yann@1250
|
206 |
+ (g, "2.66666 66666 66666 66666 66666 66666 66666 667", medium_prec);
|
yann@1250
|
207 |
}
|
yann@1250
|
208 |
{
|
yann@1250
|
209 |
mpf_class f(15.0, large_prec);
|
yann@1250
|
210 |
@@ -196,7 +198,8 @@
|
yann@1250
|
211 |
mpf_class g(0.0, very_large_prec);
|
yann@1250
|
212 |
g = mpf_class(f + q);
|
yann@1250
|
213 |
ASSERT_ALWAYS_PREC
|
yann@1250
|
214 |
- (g, "15.33333 33333 33333 33333 33333 33333 33333 33333 33333 3333",
|
yann@1250
|
215 |
+ (g, "15.33333 33333 33333 33333 33333 33333 33333 33333 33333 33333"
|
yann@1250
|
216 |
+ " 33333 33",
|
yann@1250
|
217 |
large_prec);
|
yann@1250
|
218 |
}
|
yann@1250
|
219 |
}
|