pyrox.dev
1diff -Naurd mpfr-3.1.3-a/PATCHES mpfr-3.1.3-b/PATCHES
2--- mpfr-3.1.3-a/PATCHES 2015-07-02 10:49:23.950112879 +0000
3+++ mpfr-3.1.3-b/PATCHES 2015-07-02 10:49:24.042113845 +0000
4@@ -0,0 +1 @@
5+lngamma-and-doc
6diff -Naurd mpfr-3.1.3-a/VERSION mpfr-3.1.3-b/VERSION
7--- mpfr-3.1.3-a/VERSION 2015-06-19 19:55:09.000000000 +0000
8+++ mpfr-3.1.3-b/VERSION 2015-07-02 10:49:24.042113845 +0000
9@@ -1 +1 @@
10-3.1.3
11+3.1.3-p1
12diff -Naurd mpfr-3.1.3-a/doc/mpfr.texi mpfr-3.1.3-b/doc/mpfr.texi
13--- mpfr-3.1.3-a/doc/mpfr.texi 2015-06-19 19:55:11.000000000 +0000
14+++ mpfr-3.1.3-b/doc/mpfr.texi 2015-07-02 10:49:24.018113593 +0000
15@@ -810,13 +810,17 @@
16 When the input point is in the closure of the domain of the mathematical
17 function and an input argument is +0 (resp.@: @minus{}0), one considers
18 the limit when the corresponding argument approaches 0 from above
19-(resp.@: below). If the limit is not defined (e.g., @code{mpfr_log} on
20-@minus{}0), the behavior is specified in the description of the MPFR function.
21+(resp.@: below), if possible. If the limit is not defined (e.g.,
22+@code{mpfr_sqrt} and @code{mpfr_log} on @minus{}0), the behavior is
23+specified in the description of the MPFR function, but must be consistent
24+with the rule from the above paragraph (e.g., @code{mpfr_log} on @pom{}0
25+gives @minus{}Inf).
26
27 When the result is equal to 0, its sign is determined by considering the
28 limit as if the input point were not in the domain: If one approaches 0
29 from above (resp.@: below), the result is +0 (resp.@: @minus{}0);
30-for example, @code{mpfr_sin} on +0 gives +0.
31+for example, @code{mpfr_sin} on @minus{}0 gives @minus{}0 and
32+@code{mpfr_acos} on 1 gives +0 (in all rounding modes).
33 In the other cases, the sign is specified in the description of the MPFR
34 function; for example @code{mpfr_max} on @minus{}0 and +0 gives +0.
35
36@@ -832,8 +836,8 @@
37 @c that advantages in practice), like for any bug fix.
38 Example: @code{mpfr_hypot} on (NaN,0) gives NaN, but @code{mpfr_hypot}
39 on (NaN,+Inf) gives +Inf (as specified in @ref{Special Functions}),
40-since for any finite input @var{x}, @code{mpfr_hypot} on (@var{x},+Inf)
41-gives +Inf.
42+since for any finite or infinite input @var{x}, @code{mpfr_hypot} on
43+(@var{x},+Inf) gives +Inf.
44
45 @node Exceptions, Memory Handling, Floating-Point Values on Special Numbers, MPFR Basics
46 @comment node-name, next, previous, up
47@@ -1581,7 +1585,8 @@
48 @deftypefunx int mpfr_add_z (mpfr_t @var{rop}, mpfr_t @var{op1}, mpz_t @var{op2}, mpfr_rnd_t @var{rnd})
49 @deftypefunx int mpfr_add_q (mpfr_t @var{rop}, mpfr_t @var{op1}, mpq_t @var{op2}, mpfr_rnd_t @var{rnd})
50 Set @var{rop} to @math{@var{op1} + @var{op2}} rounded in the direction
51-@var{rnd}. For types having no signed zero, it is considered unsigned
52+@var{rnd}. The IEEE-754 rules are used, in particular for signed zeros.
53+But for types having no signed zeros, 0 is considered unsigned
54 (i.e., (+0) + 0 = (+0) and (@minus{}0) + 0 = (@minus{}0)).
55 The @code{mpfr_add_d} function assumes that the radix of the @code{double} type
56 is a power of 2, with a precision at most that declared by the C implementation
57@@ -1599,7 +1604,8 @@
58 @deftypefunx int mpfr_sub_z (mpfr_t @var{rop}, mpfr_t @var{op1}, mpz_t @var{op2}, mpfr_rnd_t @var{rnd})
59 @deftypefunx int mpfr_sub_q (mpfr_t @var{rop}, mpfr_t @var{op1}, mpq_t @var{op2}, mpfr_rnd_t @var{rnd})
60 Set @var{rop} to @math{@var{op1} - @var{op2}} rounded in the direction
61-@var{rnd}. For types having no signed zero, it is considered unsigned
62+@var{rnd}. The IEEE-754 rules are used, in particular for signed zeros.
63+But for types having no signed zeros, 0 is considered unsigned
64 (i.e., (+0) @minus{} 0 = (+0), (@minus{}0) @minus{} 0 = (@minus{}0),
65 0 @minus{} (+0) = (@minus{}0) and 0 @minus{} (@minus{}0) = (+0)).
66 The same restrictions than for @code{mpfr_add_d} apply to @code{mpfr_d_sub}
67@@ -1615,7 +1621,7 @@
68 Set @var{rop} to @math{@var{op1} @GMPtimes{} @var{op2}} rounded in the
69 direction @var{rnd}.
70 When a result is zero, its sign is the product of the signs of the operands
71-(for types having no signed zero, it is considered positive).
72+(for types having no signed zeros, 0 is considered positive).
73 The same restrictions than for @code{mpfr_add_d} apply to @code{mpfr_mul_d}.
74 @end deftypefun
75
76@@ -1635,7 +1641,7 @@
77 @deftypefunx int mpfr_div_q (mpfr_t @var{rop}, mpfr_t @var{op1}, mpq_t @var{op2}, mpfr_rnd_t @var{rnd})
78 Set @var{rop} to @math{@var{op1}/@var{op2}} rounded in the direction @var{rnd}.
79 When a result is zero, its sign is the product of the signs of the operands
80-(for types having no signed zero, it is considered positive).
81+(for types having no signed zeros, 0 is considered positive).
82 The same restrictions than for @code{mpfr_add_d} apply to @code{mpfr_d_div}
83 and @code{mpfr_div_d}.
84 @end deftypefun
85@@ -1643,15 +1649,18 @@
86 @deftypefun int mpfr_sqrt (mpfr_t @var{rop}, mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
87 @deftypefunx int mpfr_sqrt_ui (mpfr_t @var{rop}, unsigned long int @var{op}, mpfr_rnd_t @var{rnd})
88 Set @var{rop} to @m{\sqrt{@var{op}}, the square root of @var{op}}
89-rounded in the direction @var{rnd} (set @var{rop} to @minus{}0 if @var{op} is
90-@minus{}0, to be consistent with the IEEE 754 standard).
91+rounded in the direction @var{rnd}. Set @var{rop} to @minus{}0 if
92+@var{op} is @minus{}0, to be consistent with the IEEE 754 standard.
93 Set @var{rop} to NaN if @var{op} is negative.
94 @end deftypefun
95
96 @deftypefun int mpfr_rec_sqrt (mpfr_t @var{rop}, mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
97 Set @var{rop} to @m{1/\sqrt{@var{op}}, the reciprocal square root of @var{op}}
98-rounded in the direction @var{rnd}. Set @var{rop} to +Inf if @var{op} is
99-@pom{}0, +0 if @var{op} is +Inf, and NaN if @var{op} is negative.
100+rounded in the direction @var{rnd}. Set @var{rop} to +Inf if @var{op} is
101+@pom{}0, +0 if @var{op} is +Inf, and NaN if @var{op} is negative. Warning!
102+Therefore the result on @minus{}0 is different from the one of the rSqrt
103+function recommended by the IEEE 754-2008 standard (Section 9.2.1), which
104+is @minus{}Inf instead of +Inf.
105 @end deftypefun
106
107 @deftypefun int mpfr_cbrt (mpfr_t @var{rop}, mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
108@@ -1832,7 +1841,9 @@
109 @m{\log_2 @var{op}, log2(@var{op})} or
110 @m{\log_{10} @var{op}, log10(@var{op})}, respectively,
111 rounded in the direction @var{rnd}.
112-Set @var{rop} to @minus{}Inf if @var{op} is @minus{}0
113+Set @var{rop} to +0 if @var{op} is 1 (in all rounding modes),
114+for consistency with the ISO C99 and IEEE 754-2008 standards.
115+Set @var{rop} to @minus{}Inf if @var{op} is @pom{}0
116 (i.e., the sign of the zero has no influence on the result).
117 @end deftypefun
118
119@@ -2003,8 +2014,11 @@
120 @deftypefun int mpfr_lngamma (mpfr_t @var{rop}, mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
121 Set @var{rop} to the value of the logarithm of the Gamma function on @var{op},
122 rounded in the direction @var{rnd}.
123-When @math{@minus{}2@var{k}@minus{}1 @le{} @var{op} @le{} @minus{}2@var{k}},
124-@var{k} being a non-negative integer, @var{rop} is set to NaN.
125+When @var{op} is 1 or 2, set @var{rop} to +0 (in all rounding modes).
126+When @var{op} is an infinity or a nonpositive integer, set @var{rop} to +Inf,
127+following the general rules on special values.
128+When @math{@minus{}2@var{k}@minus{}1 < @var{op} < @minus{}2@var{k}},
129+@var{k} being a nonnegative integer, set @var{rop} to NaN@.
130 See also @code{mpfr_lgamma}.
131 @end deftypefun
132
133@@ -2012,10 +2026,11 @@
134 Set @var{rop} to the value of the logarithm of the absolute value of the
135 Gamma function on @var{op}, rounded in the direction @var{rnd}. The sign
136 (1 or @minus{}1) of Gamma(@var{op}) is returned in the object pointed to
137-by @var{signp}. When @var{op} is an infinity or a non-positive integer, set
138-@var{rop} to +Inf. When @var{op} is NaN, @minus{}Inf or a negative integer,
139-*@var{signp} is undefined, and when @var{op} is @pom{}0, *@var{signp} is
140-the sign of the zero.
141+by @var{signp}.
142+When @var{op} is 1 or 2, set @var{rop} to +0 (in all rounding modes).
143+When @var{op} is an infinity or a nonpositive integer, set @var{rop} to +Inf.
144+When @var{op} is NaN, @minus{}Inf or a negative integer, *@var{signp} is
145+undefined, and when @var{op} is @pom{}0, *@var{signp} is the sign of the zero.
146 @end deftypefun
147
148 @deftypefun int mpfr_digamma (mpfr_t @var{rop}, mpfr_t @var{op}, mpfr_rnd_t @var{rnd})
149@@ -2064,7 +2079,10 @@
150 @deftypefunx int mpfr_fms (mpfr_t @var{rop}, mpfr_t @var{op1}, mpfr_t @var{op2}, mpfr_t @var{op3}, mpfr_rnd_t @var{rnd})
151 Set @var{rop} to @math{(@var{op1} @GMPtimes{} @var{op2}) + @var{op3}}
152 (resp.@: @math{(@var{op1} @GMPtimes{} @var{op2}) - @var{op3}})
153-rounded in the direction @var{rnd}.
154+rounded in the direction @var{rnd}. Concerning special values (signed zeros,
155+infinities, NaN), these functions behave like a multiplication followed by a
156+separate addition or subtraction. That is, the fused operation matters only
157+for rounding.
158 @end deftypefun
159
160 @deftypefun int mpfr_agm (mpfr_t @var{rop}, mpfr_t @var{op1}, mpfr_t @var{op2}, mpfr_rnd_t @var{rnd})
161@@ -2089,8 +2107,8 @@
162 i.e., $\sqrt{x^2+y^2}$,
163 @end tex
164 rounded in the direction @var{rnd}.
165-Special values are handled as described in Section F.9.4.3 of
166-the ISO C99 and IEEE 754-2008 standards:
167+Special values are handled as described in the ISO C99 (Section F.9.4.3)
168+and IEEE 754-2008 (Section 9.2.1) standards:
169 If @var{x} or @var{y} is an infinity, then +Inf is returned in @var{rop},
170 even if the other number is NaN.
171 @end deftypefun
172diff -Naurd mpfr-3.1.3-a/doc/mpfr.info mpfr-3.1.3-b/doc/mpfr.info
173--- mpfr-3.1.3-a/doc/mpfr.info 2015-06-19 19:55:53.000000000 +0000
174+++ mpfr-3.1.3-b/doc/mpfr.info 2015-07-02 10:49:38.718267817 +0000
175@@ -1,4 +1,4 @@
176-This is mpfr.info, produced by makeinfo version 5.2 from mpfr.texi.
177+This is mpfr.info, produced by makeinfo version 6.0 from mpfr.texi.
178
179 This manual documents how to install and use the Multiple Precision
180 Floating-Point Reliable Library, version 3.1.3.
181@@ -55,7 +55,7 @@
182 MPFR Copying Conditions
183 ***********************
184
185-The GNU MPFR library (or MPFR for short) is "free"; this means that
186+The GNU MPFR library (or MPFR for short) is “free”; this means that
187 everyone is free to use it and free to redistribute it on a free basis.
188 The library is not in the public domain; it is copyrighted and there are
189 restrictions on its distribution, but these restrictions are designed to
190@@ -418,7 +418,7 @@
191 4.2 Nomenclature and Types
192 ==========================
193
194-A "floating-point number", or "float" for short, is an arbitrary
195+A “floating-point number”, or “float” for short, is an arbitrary
196 precision significand (also called mantissa) with a limited precision
197 exponent. The C data type for such objects is ‘mpfr_t’ (internally
198 defined as a one-element array of a structure, and ‘mpfr_ptr’ is the C
199@@ -432,7 +432,7 @@
200 to the other functions supported by MPFR. Unless documented otherwise,
201 the sign bit of a NaN is unspecified.
202
203-The "precision" is the number of bits used to represent the significand
204+The “precision” is the number of bits used to represent the significand
205 of a floating-point number; the corresponding C data type is
206 ‘mpfr_prec_t’. The precision can be any integer between ‘MPFR_PREC_MIN’
207 and ‘MPFR_PREC_MAX’. In the current implementation, ‘MPFR_PREC_MIN’ is
208@@ -446,7 +446,7 @@
209 may abort, crash or have undefined behavior (depending on your C
210 implementation).
211
212-The "rounding mode" specifies the way to round the result of a
213+The “rounding mode” specifies the way to round the result of a
214 floating-point operation, in case the exact result can not be
215 represented exactly in the destination significand; the corresponding C
216 data type is ‘mpfr_rnd_t’.
217@@ -499,14 +499,14 @@
218 representable numbers, it is rounded to the one with the least
219 significant bit set to zero. For example, the number 2.5, which is
220 represented by (10.1) in binary, is rounded to (10.0)=2 with a precision
221-of two bits, and not to (11.0)=3. This rule avoids the "drift"
222+of two bits, and not to (11.0)=3. This rule avoids the “drift”
223 phenomenon mentioned by Knuth in volume 2 of The Art of Computer
224 Programming (Section 4.2.2).
225
226 Most MPFR functions take as first argument the destination variable,
227 as second and following arguments the input variables, as last argument
228 a rounding mode, and have a return value of type ‘int’, called the
229-"ternary value". The value stored in the destination variable is
230+“ternary value”. The value stored in the destination variable is
231 correctly rounded, i.e., MPFR behaves as if it computed the result with
232 an infinite precision, then rounded it to the precision of this
233 variable. The input variables are regarded as exact (in particular,
234@@ -572,15 +572,18 @@
235 When the input point is in the closure of the domain of the
236 mathematical function and an input argument is +0 (resp. −0), one
237 considers the limit when the corresponding argument approaches 0 from
238-above (resp. below). If the limit is not defined (e.g., ‘mpfr_log’ on
239-−0), the behavior is specified in the description of the MPFR function.
240+above (resp. below), if possible. If the limit is not defined (e.g.,
241+‘mpfr_sqrt’ and ‘mpfr_log’ on −0), the behavior is specified in the
242+description of the MPFR function, but must be consistent with the rule
243+from the above paragraph (e.g., ‘mpfr_log’ on ±0 gives −Inf).
244
245 When the result is equal to 0, its sign is determined by considering
246 the limit as if the input point were not in the domain: If one
247 approaches 0 from above (resp. below), the result is +0 (resp. −0); for
248-example, ‘mpfr_sin’ on +0 gives +0. In the other cases, the sign is
249-specified in the description of the MPFR function; for example
250-‘mpfr_max’ on −0 and +0 gives +0.
251+example, ‘mpfr_sin’ on −0 gives −0 and ‘mpfr_acos’ on 1 gives +0 (in all
252+rounding modes). In the other cases, the sign is specified in the
253+description of the MPFR function; for example ‘mpfr_max’ on −0 and +0
254+gives +0.
255
256 When the input point is not in the closure of the domain of the
257 function, the result is NaN. Example: ‘mpfr_sqrt’ on −17 gives NaN.
258@@ -590,8 +593,8 @@
259 numbers; such a case is always explicitly specified in *note MPFR
260 Interface::. Example: ‘mpfr_hypot’ on (NaN,0) gives NaN, but
261 ‘mpfr_hypot’ on (NaN,+Inf) gives +Inf (as specified in *note Special
262-Functions::), since for any finite input X, ‘mpfr_hypot’ on (X,+Inf)
263-gives +Inf.
264+Functions::), since for any finite or infinite input X, ‘mpfr_hypot’ on
265+(X,+Inf) gives +Inf.
266
267
268 File: mpfr.info, Node: Exceptions, Next: Memory Handling, Prev: Floating-Point Values on Special Numbers, Up: MPFR Basics
269@@ -1253,8 +1256,9 @@
270 mpfr_rnd_t RND)
271 -- Function: int mpfr_add_q (mpfr_t ROP, mpfr_t OP1, mpq_t OP2,
272 mpfr_rnd_t RND)
273- Set ROP to OP1 + OP2 rounded in the direction RND. For types
274- having no signed zero, it is considered unsigned (i.e., (+0) + 0 =
275+ Set ROP to OP1 + OP2 rounded in the direction RND. The IEEE-754
276+ rules are used, in particular for signed zeros. But for types
277+ having no signed zeros, 0 is considered unsigned (i.e., (+0) + 0 =
278 (+0) and (−0) + 0 = (−0)). The ‘mpfr_add_d’ function assumes that
279 the radix of the ‘double’ type is a power of 2, with a precision at
280 most that declared by the C implementation (macro
281@@ -1280,8 +1284,9 @@
282 mpfr_rnd_t RND)
283 -- Function: int mpfr_sub_q (mpfr_t ROP, mpfr_t OP1, mpq_t OP2,
284 mpfr_rnd_t RND)
285- Set ROP to OP1 - OP2 rounded in the direction RND. For types
286- having no signed zero, it is considered unsigned (i.e., (+0) − 0 =
287+ Set ROP to OP1 - OP2 rounded in the direction RND. The IEEE-754
288+ rules are used, in particular for signed zeros. But for types
289+ having no signed zeros, 0 is considered unsigned (i.e., (+0) − 0 =
290 (+0), (−0) − 0 = (−0), 0 − (+0) = (−0) and 0 − (−0) = (+0)). The
291 same restrictions than for ‘mpfr_add_d’ apply to ‘mpfr_d_sub’ and
292 ‘mpfr_sub_d’.
293@@ -1300,7 +1305,7 @@
294 mpfr_rnd_t RND)
295 Set ROP to OP1 times OP2 rounded in the direction RND. When a
296 result is zero, its sign is the product of the signs of the
297- operands (for types having no signed zero, it is considered
298+ operands (for types having no signed zeros, 0 is considered
299 positive). The same restrictions than for ‘mpfr_add_d’ apply to
300 ‘mpfr_mul_d’.
301
302@@ -1327,21 +1332,24 @@
303 mpfr_rnd_t RND)
304 Set ROP to OP1/OP2 rounded in the direction RND. When a result is
305 zero, its sign is the product of the signs of the operands (for
306- types having no signed zero, it is considered positive). The same
307+ types having no signed zeros, 0 is considered positive). The same
308 restrictions than for ‘mpfr_add_d’ apply to ‘mpfr_d_div’ and
309 ‘mpfr_div_d’.
310
311 -- Function: int mpfr_sqrt (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
312 -- Function: int mpfr_sqrt_ui (mpfr_t ROP, unsigned long int OP,
313 mpfr_rnd_t RND)
314- Set ROP to the square root of OP rounded in the direction RND (set
315- ROP to −0 if OP is −0, to be consistent with the IEEE 754
316- standard). Set ROP to NaN if OP is negative.
317+ Set ROP to the square root of OP rounded in the direction RND. Set
318+ ROP to −0 if OP is −0, to be consistent with the IEEE 754 standard.
319+ Set ROP to NaN if OP is negative.
320
321 -- Function: int mpfr_rec_sqrt (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
322 Set ROP to the reciprocal square root of OP rounded in the
323 direction RND. Set ROP to +Inf if OP is ±0, +0 if OP is +Inf, and
324- NaN if OP is negative.
325+ NaN if OP is negative. Warning! Therefore the result on −0 is
326+ different from the one of the rSqrt function recommended by the
327+ IEEE 754-2008 standard (Section 9.2.1), which is −Inf instead of
328+ +Inf.
329
330 -- Function: int mpfr_cbrt (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
331 -- Function: int mpfr_root (mpfr_t ROP, mpfr_t OP, unsigned long int K,
332@@ -1515,8 +1523,10 @@
333 -- Function: int mpfr_log2 (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
334 -- Function: int mpfr_log10 (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
335 Set ROP to the natural logarithm of OP, log2(OP) or log10(OP),
336- respectively, rounded in the direction RND. Set ROP to −Inf if OP
337- is −0 (i.e., the sign of the zero has no influence on the result).
338+ respectively, rounded in the direction RND. Set ROP to +0 if OP is
339+ 1 (in all rounding modes), for consistency with the ISO C99 and
340+ IEEE 754-2008 standards. Set ROP to −Inf if OP is ±0 (i.e., the
341+ sign of the zero has no influence on the result).
342
343 -- Function: int mpfr_exp (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
344 -- Function: int mpfr_exp2 (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
345@@ -1649,17 +1659,21 @@
346
347 -- Function: int mpfr_lngamma (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
348 Set ROP to the value of the logarithm of the Gamma function on OP,
349- rounded in the direction RND. When −2K−1 <= OP <= −2K, K being a
350- non-negative integer, ROP is set to NaN. See also ‘mpfr_lgamma’.
351+ rounded in the direction RND. When OP is 1 or 2, set ROP to +0 (in
352+ all rounding modes). When OP is an infinity or a nonpositive
353+ integer, set ROP to +Inf, following the general rules on special
354+ values. When −2K−1 < OP < −2K, K being a nonnegative integer, set
355+ ROP to NaN. See also ‘mpfr_lgamma’.
356
357 -- Function: int mpfr_lgamma (mpfr_t ROP, int *SIGNP, mpfr_t OP,
358 mpfr_rnd_t RND)
359 Set ROP to the value of the logarithm of the absolute value of the
360 Gamma function on OP, rounded in the direction RND. The sign (1 or
361 −1) of Gamma(OP) is returned in the object pointed to by SIGNP.
362- When OP is an infinity or a non-positive integer, set ROP to +Inf.
363- When OP is NaN, −Inf or a negative integer, *SIGNP is undefined,
364- and when OP is ±0, *SIGNP is the sign of the zero.
365+ When OP is 1 or 2, set ROP to +0 (in all rounding modes). When OP
366+ is an infinity or a nonpositive integer, set ROP to +Inf. When OP
367+ is NaN, −Inf or a negative integer, *SIGNP is undefined, and when
368+ OP is ±0, *SIGNP is the sign of the zero.
369
370 -- Function: int mpfr_digamma (mpfr_t ROP, mpfr_t OP, mpfr_rnd_t RND)
371 Set ROP to the value of the Digamma (sometimes also called Psi)
372@@ -1703,7 +1717,10 @@
373 -- Function: int mpfr_fms (mpfr_t ROP, mpfr_t OP1, mpfr_t OP2, mpfr_t
374 OP3, mpfr_rnd_t RND)
375 Set ROP to (OP1 times OP2) + OP3 (resp. (OP1 times OP2) - OP3)
376- rounded in the direction RND.
377+ rounded in the direction RND. Concerning special values (signed
378+ zeros, infinities, NaN), these functions behave like a
379+ multiplication followed by a separate addition or subtraction.
380+ That is, the fused operation matters only for rounding.
381
382 -- Function: int mpfr_agm (mpfr_t ROP, mpfr_t OP1, mpfr_t OP2,
383 mpfr_rnd_t RND)
384@@ -1717,9 +1734,10 @@
385 RND)
386 Set ROP to the Euclidean norm of X and Y, i.e., the square root of
387 the sum of the squares of X and Y, rounded in the direction RND.
388- Special values are handled as described in Section F.9.4.3 of the
389- ISO C99 and IEEE 754-2008 standards: If X or Y is an infinity, then
390- +Inf is returned in ROP, even if the other number is NaN.
391+ Special values are handled as described in the ISO C99 (Section
392+ F.9.4.3) and IEEE 754-2008 (Section 9.2.1) standards: If X or Y is
393+ an infinity, then +Inf is returned in ROP, even if the other number
394+ is NaN.
395
396 -- Function: int mpfr_ai (mpfr_t ROP, mpfr_t X, mpfr_rnd_t RND)
397 Set ROP to the value of the Airy function Ai on X, rounded in the
398@@ -2670,7 +2688,7 @@
399 5.16 Internals
400 ==============
401
402-A "limb" means the part of a multi-precision number that fits in a
403+A “limb” means the part of a multi-precision number that fits in a
404 single word. Usually a limb contains 32 or 64 bits. The C data type
405 for a limb is ‘mp_limb_t’.
406
407@@ -3140,7 +3158,7 @@
408 0. PREAMBLE
409
410 The purpose of this License is to make a manual, textbook, or other
411- functional and useful document "free" in the sense of freedom: to
412+ functional and useful document “free” in the sense of freedom: to
413 assure everyone the effective freedom to copy and redistribute it,
414 with or without modifying it, either commercially or
415 noncommercially. Secondarily, this License preserves for the
416@@ -3655,9 +3673,9 @@
417 * Menu:
418
419 * mpfr_abs: Basic Arithmetic Functions.
420- (line 160)
421-* mpfr_acos: Special Functions. (line 51)
422-* mpfr_acosh: Special Functions. (line 115)
423+ (line 165)
424+* mpfr_acos: Special Functions. (line 53)
425+* mpfr_acosh: Special Functions. (line 117)
426 * mpfr_add: Basic Arithmetic Functions.
427 (line 6)
428 * mpfr_add_d: Basic Arithmetic Functions.
429@@ -3670,15 +3688,15 @@
430 (line 8)
431 * mpfr_add_z: Basic Arithmetic Functions.
432 (line 14)
433-* mpfr_agm: Special Functions. (line 210)
434-* mpfr_ai: Special Functions. (line 226)
435-* mpfr_asin: Special Functions. (line 52)
436-* mpfr_asinh: Special Functions. (line 116)
437+* mpfr_agm: Special Functions. (line 219)
438+* mpfr_ai: Special Functions. (line 236)
439+* mpfr_asin: Special Functions. (line 54)
440+* mpfr_asinh: Special Functions. (line 118)
441 * mpfr_asprintf: Formatted Output Functions.
442 (line 193)
443-* mpfr_atan: Special Functions. (line 53)
444-* mpfr_atan2: Special Functions. (line 63)
445-* mpfr_atanh: Special Functions. (line 117)
446+* mpfr_atan: Special Functions. (line 55)
447+* mpfr_atan2: Special Functions. (line 65)
448+* mpfr_atanh: Special Functions. (line 119)
449 * mpfr_buildopt_decimal_p: Miscellaneous Functions.
450 (line 162)
451 * mpfr_buildopt_gmpinternals_p: Miscellaneous Functions.
452@@ -3690,7 +3708,7 @@
453 * mpfr_can_round: Rounding Related Functions.
454 (line 39)
455 * mpfr_cbrt: Basic Arithmetic Functions.
456- (line 108)
457+ (line 113)
458 * mpfr_ceil: Integer Related Functions.
459 (line 7)
460 * mpfr_check_range: Exception Related Functions.
461@@ -3735,18 +3753,18 @@
462 (line 27)
463 * mpfr_cmp_z: Comparison Functions.
464 (line 11)
465-* mpfr_const_catalan: Special Functions. (line 237)
466-* mpfr_const_euler: Special Functions. (line 236)
467-* mpfr_const_log2: Special Functions. (line 234)
468-* mpfr_const_pi: Special Functions. (line 235)
469+* mpfr_const_catalan: Special Functions. (line 247)
470+* mpfr_const_euler: Special Functions. (line 246)
471+* mpfr_const_log2: Special Functions. (line 244)
472+* mpfr_const_pi: Special Functions. (line 245)
473 * mpfr_copysign: Miscellaneous Functions.
474 (line 109)
475-* mpfr_cos: Special Functions. (line 29)
476-* mpfr_cosh: Special Functions. (line 95)
477-* mpfr_cot: Special Functions. (line 47)
478-* mpfr_coth: Special Functions. (line 111)
479-* mpfr_csc: Special Functions. (line 46)
480-* mpfr_csch: Special Functions. (line 110)
481+* mpfr_cos: Special Functions. (line 31)
482+* mpfr_cosh: Special Functions. (line 97)
483+* mpfr_cot: Special Functions. (line 49)
484+* mpfr_coth: Special Functions. (line 113)
485+* mpfr_csc: Special Functions. (line 48)
486+* mpfr_csch: Special Functions. (line 112)
487 * mpfr_custom_get_exp: Custom Interface. (line 75)
488 * mpfr_custom_get_kind: Custom Interface. (line 65)
489 * mpfr_custom_get_significand: Custom Interface. (line 70)
490@@ -3756,47 +3774,47 @@
491 * mpfr_custom_move: Custom Interface. (line 82)
492 * MPFR_DECL_INIT: Initialization Functions.
493 (line 74)
494-* mpfr_digamma: Special Functions. (line 166)
495+* mpfr_digamma: Special Functions. (line 172)
496 * mpfr_dim: Basic Arithmetic Functions.
497- (line 166)
498+ (line 171)
499 * mpfr_div: Basic Arithmetic Functions.
500- (line 72)
501+ (line 74)
502 * mpfr_divby0_p: Exception Related Functions.
503 (line 134)
504 * mpfr_div_2exp: Compatibility with MPF.
505 (line 49)
506 * mpfr_div_2si: Basic Arithmetic Functions.
507- (line 181)
508+ (line 186)
509 * mpfr_div_2ui: Basic Arithmetic Functions.
510- (line 179)
511+ (line 184)
512 * mpfr_div_d: Basic Arithmetic Functions.
513- (line 84)
514+ (line 86)
515 * mpfr_div_q: Basic Arithmetic Functions.
516- (line 88)
517+ (line 90)
518 * mpfr_div_si: Basic Arithmetic Functions.
519- (line 80)
520+ (line 82)
521 * mpfr_div_ui: Basic Arithmetic Functions.
522- (line 76)
523+ (line 78)
524 * mpfr_div_z: Basic Arithmetic Functions.
525- (line 86)
526+ (line 88)
527 * mpfr_d_div: Basic Arithmetic Functions.
528- (line 82)
529+ (line 84)
530 * mpfr_d_sub: Basic Arithmetic Functions.
531- (line 35)
532-* mpfr_eint: Special Functions. (line 133)
533+ (line 36)
534+* mpfr_eint: Special Functions. (line 135)
535 * mpfr_eq: Compatibility with MPF.
536 (line 28)
537 * mpfr_equal_p: Comparison Functions.
538 (line 59)
539 * mpfr_erangeflag_p: Exception Related Functions.
540 (line 137)
541-* mpfr_erf: Special Functions. (line 177)
542-* mpfr_erfc: Special Functions. (line 178)
543-* mpfr_exp: Special Functions. (line 23)
544-* mpfr_exp10: Special Functions. (line 25)
545-* mpfr_exp2: Special Functions. (line 24)
546-* mpfr_expm1: Special Functions. (line 129)
547-* mpfr_fac_ui: Special Functions. (line 121)
548+* mpfr_erf: Special Functions. (line 183)
549+* mpfr_erfc: Special Functions. (line 184)
550+* mpfr_exp: Special Functions. (line 25)
551+* mpfr_exp10: Special Functions. (line 27)
552+* mpfr_exp2: Special Functions. (line 26)
553+* mpfr_expm1: Special Functions. (line 131)
554+* mpfr_fac_ui: Special Functions. (line 123)
555 * mpfr_fits_intmax_p: Conversion Functions.
556 (line 150)
557 * mpfr_fits_sint_p: Conversion Functions.
558@@ -3815,20 +3833,20 @@
559 (line 147)
560 * mpfr_floor: Integer Related Functions.
561 (line 8)
562-* mpfr_fma: Special Functions. (line 203)
563+* mpfr_fma: Special Functions. (line 209)
564 * mpfr_fmod: Integer Related Functions.
565 (line 92)
566-* mpfr_fms: Special Functions. (line 205)
567+* mpfr_fms: Special Functions. (line 211)
568 * mpfr_fprintf: Formatted Output Functions.
569 (line 157)
570 * mpfr_frac: Integer Related Functions.
571 (line 76)
572-* mpfr_free_cache: Special Functions. (line 244)
573+* mpfr_free_cache: Special Functions. (line 254)
574 * mpfr_free_str: Conversion Functions.
575 (line 137)
576 * mpfr_frexp: Conversion Functions.
577 (line 45)
578-* mpfr_gamma: Special Functions. (line 148)
579+* mpfr_gamma: Special Functions. (line 150)
580 * mpfr_get_d: Conversion Functions.
581 (line 7)
582 * mpfr_get_decimal64: Conversion Functions.
583@@ -3887,7 +3905,7 @@
584 (line 56)
585 * mpfr_greater_p: Comparison Functions.
586 (line 55)
587-* mpfr_hypot: Special Functions. (line 218)
588+* mpfr_hypot: Special Functions. (line 227)
589 * mpfr_inexflag_p: Exception Related Functions.
590 (line 136)
591 * mpfr_inf_p: Comparison Functions.
592@@ -3922,21 +3940,21 @@
593 (line 31)
594 * mpfr_integer_p: Integer Related Functions.
595 (line 119)
596-* mpfr_j0: Special Functions. (line 182)
597-* mpfr_j1: Special Functions. (line 183)
598-* mpfr_jn: Special Functions. (line 184)
599+* mpfr_j0: Special Functions. (line 188)
600+* mpfr_j1: Special Functions. (line 189)
601+* mpfr_jn: Special Functions. (line 190)
602 * mpfr_lessequal_p: Comparison Functions.
603 (line 58)
604 * mpfr_lessgreater_p: Comparison Functions.
605 (line 64)
606 * mpfr_less_p: Comparison Functions.
607 (line 57)
608-* mpfr_lgamma: Special Functions. (line 157)
609-* mpfr_li2: Special Functions. (line 143)
610-* mpfr_lngamma: Special Functions. (line 152)
611+* mpfr_lgamma: Special Functions. (line 162)
612+* mpfr_li2: Special Functions. (line 145)
613+* mpfr_lngamma: Special Functions. (line 154)
614 * mpfr_log: Special Functions. (line 16)
615 * mpfr_log10: Special Functions. (line 18)
616-* mpfr_log1p: Special Functions. (line 125)
617+* mpfr_log1p: Special Functions. (line 127)
618 * mpfr_log2: Special Functions. (line 17)
619 * mpfr_max: Miscellaneous Functions.
620 (line 22)
621@@ -3947,29 +3965,29 @@
622 * mpfr_modf: Integer Related Functions.
623 (line 82)
624 * mpfr_mul: Basic Arithmetic Functions.
625- (line 51)
626+ (line 53)
627 * mpfr_mul_2exp: Compatibility with MPF.
628 (line 47)
629 * mpfr_mul_2si: Basic Arithmetic Functions.
630- (line 174)
631+ (line 179)
632 * mpfr_mul_2ui: Basic Arithmetic Functions.
633- (line 172)
634+ (line 177)
635 * mpfr_mul_d: Basic Arithmetic Functions.
636- (line 57)
637+ (line 59)
638 * mpfr_mul_q: Basic Arithmetic Functions.
639- (line 61)
640+ (line 63)
641 * mpfr_mul_si: Basic Arithmetic Functions.
642- (line 55)
643+ (line 57)
644 * mpfr_mul_ui: Basic Arithmetic Functions.
645- (line 53)
646+ (line 55)
647 * mpfr_mul_z: Basic Arithmetic Functions.
648- (line 59)
649+ (line 61)
650 * mpfr_nanflag_p: Exception Related Functions.
651 (line 135)
652 * mpfr_nan_p: Comparison Functions.
653 (line 39)
654 * mpfr_neg: Basic Arithmetic Functions.
655- (line 159)
656+ (line 164)
657 * mpfr_nextabove: Miscellaneous Functions.
658 (line 15)
659 * mpfr_nextbelow: Miscellaneous Functions.
660@@ -3983,13 +4001,13 @@
661 * mpfr_overflow_p: Exception Related Functions.
662 (line 133)
663 * mpfr_pow: Basic Arithmetic Functions.
664- (line 116)
665+ (line 121)
666 * mpfr_pow_si: Basic Arithmetic Functions.
667- (line 120)
668+ (line 125)
669 * mpfr_pow_ui: Basic Arithmetic Functions.
670- (line 118)
671+ (line 123)
672 * mpfr_pow_z: Basic Arithmetic Functions.
673- (line 122)
674+ (line 127)
675 * mpfr_prec_round: Rounding Related Functions.
676 (line 13)
677 * ‘mpfr_prec_t’: Nomenclature and Types.
678@@ -3999,7 +4017,7 @@
679 * mpfr_print_rnd_mode: Rounding Related Functions.
680 (line 71)
681 * mpfr_rec_sqrt: Basic Arithmetic Functions.
682- (line 103)
683+ (line 105)
684 * mpfr_regular_p: Comparison Functions.
685 (line 43)
686 * mpfr_reldiff: Compatibility with MPF.
687@@ -4021,11 +4039,11 @@
688 * ‘mpfr_rnd_t’: Nomenclature and Types.
689 (line 34)
690 * mpfr_root: Basic Arithmetic Functions.
691- (line 109)
692+ (line 114)
693 * mpfr_round: Integer Related Functions.
694 (line 9)
695-* mpfr_sec: Special Functions. (line 45)
696-* mpfr_sech: Special Functions. (line 109)
697+* mpfr_sec: Special Functions. (line 47)
698+* mpfr_sech: Special Functions. (line 111)
699 * mpfr_set: Assignment Functions.
700 (line 9)
701 * mpfr_setsign: Miscellaneous Functions.
702@@ -4100,57 +4118,57 @@
703 (line 49)
704 * mpfr_signbit: Miscellaneous Functions.
705 (line 99)
706-* mpfr_sin: Special Functions. (line 30)
707-* mpfr_sinh: Special Functions. (line 96)
708-* mpfr_sinh_cosh: Special Functions. (line 101)
709-* mpfr_sin_cos: Special Functions. (line 35)
710+* mpfr_sin: Special Functions. (line 32)
711+* mpfr_sinh: Special Functions. (line 98)
712+* mpfr_sinh_cosh: Special Functions. (line 103)
713+* mpfr_sin_cos: Special Functions. (line 37)
714 * mpfr_si_div: Basic Arithmetic Functions.
715- (line 78)
716+ (line 80)
717 * mpfr_si_sub: Basic Arithmetic Functions.
718- (line 31)
719+ (line 32)
720 * mpfr_snprintf: Formatted Output Functions.
721 (line 180)
722 * mpfr_sprintf: Formatted Output Functions.
723 (line 170)
724 * mpfr_sqr: Basic Arithmetic Functions.
725- (line 69)
726+ (line 71)
727 * mpfr_sqrt: Basic Arithmetic Functions.
728- (line 96)
729+ (line 98)
730 * mpfr_sqrt_ui: Basic Arithmetic Functions.
731- (line 97)
732+ (line 99)
733 * mpfr_strtofr: Assignment Functions.
734 (line 80)
735 * mpfr_sub: Basic Arithmetic Functions.
736- (line 25)
737+ (line 26)
738 * mpfr_subnormalize: Exception Related Functions.
739 (line 60)
740 * mpfr_sub_d: Basic Arithmetic Functions.
741- (line 37)
742+ (line 38)
743 * mpfr_sub_q: Basic Arithmetic Functions.
744- (line 43)
745+ (line 44)
746 * mpfr_sub_si: Basic Arithmetic Functions.
747- (line 33)
748+ (line 34)
749 * mpfr_sub_ui: Basic Arithmetic Functions.
750- (line 29)
751+ (line 30)
752 * mpfr_sub_z: Basic Arithmetic Functions.
753- (line 41)
754-* mpfr_sum: Special Functions. (line 252)
755+ (line 42)
756+* mpfr_sum: Special Functions. (line 262)
757 * mpfr_swap: Assignment Functions.
758 (line 150)
759 * ‘mpfr_t’: Nomenclature and Types.
760 (line 6)
761-* mpfr_tan: Special Functions. (line 31)
762-* mpfr_tanh: Special Functions. (line 97)
763+* mpfr_tan: Special Functions. (line 33)
764+* mpfr_tanh: Special Functions. (line 99)
765 * mpfr_trunc: Integer Related Functions.
766 (line 10)
767 * mpfr_ui_div: Basic Arithmetic Functions.
768- (line 74)
769+ (line 76)
770 * mpfr_ui_pow: Basic Arithmetic Functions.
771- (line 126)
772+ (line 131)
773 * mpfr_ui_pow_ui: Basic Arithmetic Functions.
774- (line 124)
775+ (line 129)
776 * mpfr_ui_sub: Basic Arithmetic Functions.
777- (line 27)
778+ (line 28)
779 * mpfr_underflow_p: Exception Related Functions.
780 (line 132)
781 * mpfr_unordered_p: Comparison Functions.
782@@ -4181,61 +4199,61 @@
783 (line 182)
784 * mpfr_vsprintf: Formatted Output Functions.
785 (line 171)
786-* mpfr_y0: Special Functions. (line 193)
787-* mpfr_y1: Special Functions. (line 194)
788-* mpfr_yn: Special Functions. (line 195)
789+* mpfr_y0: Special Functions. (line 199)
790+* mpfr_y1: Special Functions. (line 200)
791+* mpfr_yn: Special Functions. (line 201)
792 * mpfr_zero_p: Comparison Functions.
793 (line 42)
794-* mpfr_zeta: Special Functions. (line 171)
795-* mpfr_zeta_ui: Special Functions. (line 172)
796+* mpfr_zeta: Special Functions. (line 177)
797+* mpfr_zeta_ui: Special Functions. (line 178)
798 * mpfr_z_sub: Basic Arithmetic Functions.
799- (line 39)
800+ (line 40)
801
802
803
804 Tag Table:
805 Node: Top775
806 Node: Copying2007
807-Node: Introduction to MPFR3766
808-Node: Installing MPFR5880
809-Node: Reporting Bugs11323
810-Node: MPFR Basics13353
811-Node: Headers and Libraries13669
812-Node: Nomenclature and Types16828
813-Node: MPFR Variable Conventions18874
814-Node: Rounding Modes20418
815-Ref: ternary value21544
816-Node: Floating-Point Values on Special Numbers23526
817-Node: Exceptions26572
818-Node: Memory Handling29749
819-Node: MPFR Interface30894
820-Node: Initialization Functions33008
821-Node: Assignment Functions40318
822-Node: Combined Initialization and Assignment Functions49673
823-Node: Conversion Functions50974
824-Node: Basic Arithmetic Functions60035
825-Node: Comparison Functions69200
826-Node: Special Functions72687
827-Node: Input and Output Functions86672
828-Node: Formatted Output Functions88644
829-Node: Integer Related Functions98431
830-Node: Rounding Related Functions105051
831-Node: Miscellaneous Functions108888
832-Node: Exception Related Functions117568
833-Node: Compatibility with MPF124386
834-Node: Custom Interface127127
835-Node: Internals131526
836-Node: API Compatibility133066
837-Node: Type and Macro Changes134995
838-Node: Added Functions137844
839-Node: Changed Functions141132
840-Node: Removed Functions145545
841-Node: Other Changes145973
842-Node: Contributors147576
843-Node: References150219
844-Node: GNU Free Documentation License151973
845-Node: Concept Index174562
846-Node: Function and Type Index180659
847+Node: Introduction to MPFR3770
848+Node: Installing MPFR5884
849+Node: Reporting Bugs11327
850+Node: MPFR Basics13357
851+Node: Headers and Libraries13673
852+Node: Nomenclature and Types16832
853+Node: MPFR Variable Conventions18894
854+Node: Rounding Modes20438
855+Ref: ternary value21568
856+Node: Floating-Point Values on Special Numbers23554
857+Node: Exceptions26813
858+Node: Memory Handling29990
859+Node: MPFR Interface31135
860+Node: Initialization Functions33249
861+Node: Assignment Functions40559
862+Node: Combined Initialization and Assignment Functions49914
863+Node: Conversion Functions51215
864+Node: Basic Arithmetic Functions60276
865+Node: Comparison Functions69777
866+Node: Special Functions73264
867+Node: Input and Output Functions87862
868+Node: Formatted Output Functions89834
869+Node: Integer Related Functions99621
870+Node: Rounding Related Functions106241
871+Node: Miscellaneous Functions110078
872+Node: Exception Related Functions118758
873+Node: Compatibility with MPF125576
874+Node: Custom Interface128317
875+Node: Internals132716
876+Node: API Compatibility134260
877+Node: Type and Macro Changes136189
878+Node: Added Functions139038
879+Node: Changed Functions142326
880+Node: Removed Functions146739
881+Node: Other Changes147167
882+Node: Contributors148770
883+Node: References151413
884+Node: GNU Free Documentation License153167
885+Node: Concept Index175760
886+Node: Function and Type Index181857
887
888 End Tag Table
889
890diff -Naurd mpfr-3.1.3-a/src/lngamma.c mpfr-3.1.3-b/src/lngamma.c
891--- mpfr-3.1.3-a/src/lngamma.c 2015-06-19 19:55:10.000000000 +0000
892+++ mpfr-3.1.3-b/src/lngamma.c 2015-07-02 10:49:24.018113593 +0000
893@@ -603,16 +603,17 @@
894 mpfr_get_prec (y), mpfr_log_prec, y, inex));
895
896 /* special cases */
897- if (MPFR_UNLIKELY (MPFR_IS_SINGULAR (x)))
898+ if (MPFR_UNLIKELY (MPFR_IS_SINGULAR (x) ||
899+ (MPFR_IS_NEG (x) && mpfr_integer_p (x))))
900 {
901- if (MPFR_IS_NAN (x) || MPFR_IS_NEG (x))
902+ if (MPFR_IS_NAN (x))
903 {
904 MPFR_SET_NAN (y);
905 MPFR_RET_NAN;
906 }
907- else /* lngamma(+Inf) = lngamma(+0) = +Inf */
908+ else /* lngamma(+/-Inf) = lngamma(nonpositive integer) = +Inf */
909 {
910- if (MPFR_IS_ZERO (x))
911+ if (!MPFR_IS_INF (x))
912 mpfr_set_divby0 ();
913 MPFR_SET_INF (y);
914 MPFR_SET_POS (y);
915@@ -620,8 +621,8 @@
916 }
917 }
918
919- /* if x < 0 and -2k-1 <= x <= -2k, then lngamma(x) = NaN */
920- if (MPFR_IS_NEG (x) && (unit_bit (x) == 0 || mpfr_integer_p (x)))
921+ /* if -2k-1 < x < -2k <= 0, then lngamma(x) = NaN */
922+ if (MPFR_IS_NEG (x) && unit_bit (x) == 0)
923 {
924 MPFR_SET_NAN (y);
925 MPFR_RET_NAN;
926diff -Naurd mpfr-3.1.3-a/src/mpfr.h mpfr-3.1.3-b/src/mpfr.h
927--- mpfr-3.1.3-a/src/mpfr.h 2015-06-19 19:55:10.000000000 +0000
928+++ mpfr-3.1.3-b/src/mpfr.h 2015-07-02 10:49:24.038113803 +0000
929@@ -27,7 +27,7 @@
930 #define MPFR_VERSION_MAJOR 3
931 #define MPFR_VERSION_MINOR 1
932 #define MPFR_VERSION_PATCHLEVEL 3
933-#define MPFR_VERSION_STRING "3.1.3"
934+#define MPFR_VERSION_STRING "3.1.3-p1"
935
936 /* Macros dealing with MPFR VERSION */
937 #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c))
938diff -Naurd mpfr-3.1.3-a/src/version.c mpfr-3.1.3-b/src/version.c
939--- mpfr-3.1.3-a/src/version.c 2015-06-19 19:55:10.000000000 +0000
940+++ mpfr-3.1.3-b/src/version.c 2015-07-02 10:49:24.042113845 +0000
941@@ -25,5 +25,5 @@
942 const char *
943 mpfr_get_version (void)
944 {
945- return "3.1.3";
946+ return "3.1.3-p1";
947 }
948diff -Naurd mpfr-3.1.3-a/tests/tlngamma.c mpfr-3.1.3-b/tests/tlngamma.c
949--- mpfr-3.1.3-a/tests/tlngamma.c 2015-06-19 19:55:10.000000000 +0000
950+++ mpfr-3.1.3-b/tests/tlngamma.c 2015-07-02 10:49:24.018113593 +0000
951@@ -33,7 +33,7 @@
952 special (void)
953 {
954 mpfr_t x, y;
955- int inex;
956+ int i, inex;
957
958 mpfr_init (x);
959 mpfr_init (y);
960@@ -46,25 +46,29 @@
961 exit (1);
962 }
963
964- mpfr_set_inf (x, -1);
965+ mpfr_set_inf (x, 1);
966+ mpfr_clear_flags ();
967 mpfr_lngamma (y, x, MPFR_RNDN);
968- if (!mpfr_nan_p (y))
969+ if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0 || __gmpfr_flags != 0)
970 {
971- printf ("Error for lngamma(-Inf)\n");
972+ printf ("Error for lngamma(+Inf)\n");
973 exit (1);
974 }
975
976- mpfr_set_inf (x, 1);
977+ mpfr_set_inf (x, -1);
978+ mpfr_clear_flags ();
979 mpfr_lngamma (y, x, MPFR_RNDN);
980- if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0)
981+ if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0 || __gmpfr_flags != 0)
982 {
983- printf ("Error for lngamma(+Inf)\n");
984+ printf ("Error for lngamma(-Inf)\n");
985 exit (1);
986 }
987
988 mpfr_set_ui (x, 0, MPFR_RNDN);
989+ mpfr_clear_flags ();
990 mpfr_lngamma (y, x, MPFR_RNDN);
991- if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0)
992+ if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0 ||
993+ __gmpfr_flags != MPFR_FLAGS_DIVBY0)
994 {
995 printf ("Error for lngamma(+0)\n");
996 exit (1);
997@@ -72,32 +76,58 @@
998
999 mpfr_set_ui (x, 0, MPFR_RNDN);
1000 mpfr_neg (x, x, MPFR_RNDN);
1001+ mpfr_clear_flags ();
1002 mpfr_lngamma (y, x, MPFR_RNDN);
1003- if (!mpfr_nan_p (y))
1004+ if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0 ||
1005+ __gmpfr_flags != MPFR_FLAGS_DIVBY0)
1006 {
1007 printf ("Error for lngamma(-0)\n");
1008 exit (1);
1009 }
1010
1011 mpfr_set_ui (x, 1, MPFR_RNDN);
1012+ mpfr_clear_flags ();
1013 mpfr_lngamma (y, x, MPFR_RNDN);
1014- if (MPFR_IS_NAN (y) || mpfr_cmp_ui (y, 0) || MPFR_IS_NEG (y))
1015+ if (mpfr_cmp_ui0 (y, 0) || MPFR_IS_NEG (y))
1016 {
1017 printf ("Error for lngamma(1)\n");
1018 exit (1);
1019 }
1020
1021- mpfr_set_si (x, -1, MPFR_RNDN);
1022- mpfr_lngamma (y, x, MPFR_RNDN);
1023- if (!mpfr_nan_p (y))
1024+ for (i = 1; i <= 5; i++)
1025 {
1026- printf ("Error for lngamma(-1)\n");
1027- exit (1);
1028+ int c;
1029+
1030+ mpfr_set_si (x, -i, MPFR_RNDN);
1031+ mpfr_clear_flags ();
1032+ mpfr_lngamma (y, x, MPFR_RNDN);
1033+ if (!mpfr_inf_p (y) || mpfr_sgn (y) < 0 ||
1034+ __gmpfr_flags != MPFR_FLAGS_DIVBY0)
1035+ {
1036+ printf ("Error for lngamma(-%d)\n", i);
1037+ exit (1);
1038+ }
1039+ if (i & 1)
1040+ {
1041+ mpfr_nextabove (x);
1042+ c = '+';
1043+ }
1044+ else
1045+ {
1046+ mpfr_nextbelow (x);
1047+ c = '-';
1048+ }
1049+ mpfr_lngamma (y, x, MPFR_RNDN);
1050+ if (!mpfr_nan_p (y))
1051+ {
1052+ printf ("Error for lngamma(-%d%cepsilon)\n", i, c);
1053+ exit (1);
1054+ }
1055 }
1056
1057 mpfr_set_ui (x, 2, MPFR_RNDN);
1058 mpfr_lngamma (y, x, MPFR_RNDN);
1059- if (MPFR_IS_NAN (y) || mpfr_cmp_ui (y, 0) || MPFR_IS_NEG (y))
1060+ if (mpfr_cmp_ui0 (y, 0) || MPFR_IS_NEG (y))
1061 {
1062 printf ("Error for lngamma(2)\n");
1063 exit (1);
1064@@ -127,7 +157,7 @@
1065 mpfr_set_str (x, CHECK_X2, 10, MPFR_RNDN);
1066 mpfr_lngamma (y, x, MPFR_RNDN);
1067 mpfr_set_str (x, CHECK_Y2, 10, MPFR_RNDN);
1068- if (MPFR_IS_NAN (y) || mpfr_cmp (y, x))
1069+ if (mpfr_cmp0 (y, x))
1070 {
1071 printf ("mpfr_lngamma("CHECK_X2") is wrong:\n"
1072 "expected ");
1073@@ -143,7 +173,7 @@
1074 mpfr_lngamma (y, x, MPFR_RNDU);
1075 mpfr_set_prec (x, 175);
1076 mpfr_set_str_binary (x, "0.1010001100011101101011001101110010100001000001000001110011000001101100001111001001000101011011100100010101011110100111110101010100010011010010000101010111001100011000101111E7");
1077- if (MPFR_IS_NAN (y) || mpfr_cmp (x, y))
1078+ if (mpfr_cmp0 (x, y))
1079 {
1080 printf ("Error in mpfr_lngamma (1)\n");
1081 exit (1);
1082@@ -155,7 +185,7 @@
1083 mpfr_lngamma (x, y, MPFR_RNDZ);
1084 mpfr_set_prec (y, 21);
1085 mpfr_set_str_binary (y, "0.111000101000001100101E9");
1086- if (MPFR_IS_NAN (x) || mpfr_cmp (x, y))
1087+ if (mpfr_cmp0 (x, y))
1088 {
1089 printf ("Error in mpfr_lngamma (120)\n");
1090 printf ("Expected "); mpfr_print_binary (y); puts ("");
1091@@ -169,7 +199,7 @@
1092 inex = mpfr_lngamma (y, x, MPFR_RNDN);
1093 mpfr_set_prec (x, 206);
1094 mpfr_set_str_binary (x, "0.10000111011000000011100010101001100110001110000111100011000100100110110010001011011110101001111011110110000001010100111011010000000011100110110101100111000111010011110010000100010111101010001101000110101001E13");
1095- if (MPFR_IS_NAN (y) || mpfr_cmp (x, y))
1096+ if (mpfr_cmp0 (x, y))
1097 {
1098 printf ("Error in mpfr_lngamma (768)\n");
1099 exit (1);
1100@@ -185,7 +215,7 @@
1101 mpfr_set_str_binary (x, "0.1100E-66");
1102 mpfr_lngamma (y, x, MPFR_RNDN);
1103 mpfr_set_str_binary (x, "0.1100E6");
1104- if (MPFR_IS_NAN (y) || mpfr_cmp (x, y))
1105+ if (mpfr_cmp0 (x, y))
1106 {
1107 printf ("Error for lngamma(0.1100E-66)\n");
1108 exit (1);
1109@@ -199,7 +229,7 @@
1110 mpfr_lngamma (y, x, MPFR_RNDN);
1111 mpfr_set_prec (x, 32);
1112 mpfr_set_str_binary (x, "-0.10001000111011111011000010100010E207");
1113- if (MPFR_IS_NAN (y) || mpfr_cmp (x, y))
1114+ if (mpfr_cmp0 (x, y))
1115 {
1116 printf ("Error for lngamma(-2^199+0.5)\n");
1117 printf ("Got ");
1118diff -Naurd mpfr-3.1.3-a/PATCHES mpfr-3.1.3-b/PATCHES
1119--- mpfr-3.1.3-a/PATCHES 2015-07-02 10:50:08.046573308 +0000
1120+++ mpfr-3.1.3-b/PATCHES 2015-07-02 10:50:08.126574142 +0000
1121@@ -0,0 +1 @@
1122+muldiv-2exp-overflow
1123diff -Naurd mpfr-3.1.3-a/VERSION mpfr-3.1.3-b/VERSION
1124--- mpfr-3.1.3-a/VERSION 2015-07-02 10:49:24.042113845 +0000
1125+++ mpfr-3.1.3-b/VERSION 2015-07-02 10:50:08.126574142 +0000
1126@@ -1 +1 @@
1127-3.1.3-p1
1128+3.1.3-p2
1129diff -Naurd mpfr-3.1.3-a/src/div_2si.c mpfr-3.1.3-b/src/div_2si.c
1130--- mpfr-3.1.3-a/src/div_2si.c 2015-06-19 19:55:10.000000000 +0000
1131+++ mpfr-3.1.3-b/src/div_2si.c 2015-07-02 10:50:08.106573933 +0000
1132@@ -49,7 +49,7 @@
1133 rnd_mode = MPFR_RNDZ;
1134 return mpfr_underflow (y, rnd_mode, MPFR_SIGN(y));
1135 }
1136- else if (MPFR_UNLIKELY(n < 0 && (__gmpfr_emax < MPFR_EMIN_MIN - n ||
1137+ else if (MPFR_UNLIKELY(n <= 0 && (__gmpfr_emax < MPFR_EMIN_MIN - n ||
1138 exp > __gmpfr_emax + n)) )
1139 return mpfr_overflow (y, rnd_mode, MPFR_SIGN(y));
1140
1141diff -Naurd mpfr-3.1.3-a/src/div_2ui.c mpfr-3.1.3-b/src/div_2ui.c
1142--- mpfr-3.1.3-a/src/div_2ui.c 2015-06-19 19:55:10.000000000 +0000
1143+++ mpfr-3.1.3-b/src/div_2ui.c 2015-07-02 10:50:08.106573933 +0000
1144@@ -32,7 +32,7 @@
1145 rnd_mode),
1146 ("y[%Pu]=%.*Rg inexact=%d", mpfr_get_prec(y), mpfr_log_prec, y, inexact));
1147
1148- if (MPFR_UNLIKELY (MPFR_IS_SINGULAR (x)))
1149+ if (MPFR_UNLIKELY (n == 0 || MPFR_IS_SINGULAR (x)))
1150 return mpfr_set (y, x, rnd_mode);
1151 else
1152 {
1153diff -Naurd mpfr-3.1.3-a/src/mpfr.h mpfr-3.1.3-b/src/mpfr.h
1154--- mpfr-3.1.3-a/src/mpfr.h 2015-07-02 10:49:24.038113803 +0000
1155+++ mpfr-3.1.3-b/src/mpfr.h 2015-07-02 10:50:08.126574142 +0000
1156@@ -27,7 +27,7 @@
1157 #define MPFR_VERSION_MAJOR 3
1158 #define MPFR_VERSION_MINOR 1
1159 #define MPFR_VERSION_PATCHLEVEL 3
1160-#define MPFR_VERSION_STRING "3.1.3-p1"
1161+#define MPFR_VERSION_STRING "3.1.3-p2"
1162
1163 /* Macros dealing with MPFR VERSION */
1164 #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c))
1165diff -Naurd mpfr-3.1.3-a/src/mul_2si.c mpfr-3.1.3-b/src/mul_2si.c
1166--- mpfr-3.1.3-a/src/mul_2si.c 2015-06-19 19:55:10.000000000 +0000
1167+++ mpfr-3.1.3-b/src/mul_2si.c 2015-07-02 10:50:08.106573933 +0000
1168@@ -39,7 +39,7 @@
1169 {
1170 mpfr_exp_t exp = MPFR_GET_EXP (x);
1171 MPFR_SETRAW (inexact, y, x, exp, rnd_mode);
1172- if (MPFR_UNLIKELY( n > 0 && (__gmpfr_emax < MPFR_EMIN_MIN + n ||
1173+ if (MPFR_UNLIKELY(n >= 0 && (__gmpfr_emax < MPFR_EMIN_MIN + n ||
1174 exp > __gmpfr_emax - n)))
1175 return mpfr_overflow (y, rnd_mode, MPFR_SIGN(y));
1176 else if (MPFR_UNLIKELY(n < 0 && (__gmpfr_emin > MPFR_EMAX_MAX + n ||
1177diff -Naurd mpfr-3.1.3-a/src/version.c mpfr-3.1.3-b/src/version.c
1178--- mpfr-3.1.3-a/src/version.c 2015-07-02 10:49:24.042113845 +0000
1179+++ mpfr-3.1.3-b/src/version.c 2015-07-02 10:50:08.126574142 +0000
1180@@ -25,5 +25,5 @@
1181 const char *
1182 mpfr_get_version (void)
1183 {
1184- return "3.1.3-p1";
1185+ return "3.1.3-p2";
1186 }
1187diff -Naurd mpfr-3.1.3-a/tests/tmul_2exp.c mpfr-3.1.3-b/tests/tmul_2exp.c
1188--- mpfr-3.1.3-a/tests/tmul_2exp.c 2015-06-19 19:55:10.000000000 +0000
1189+++ mpfr-3.1.3-b/tests/tmul_2exp.c 2015-07-02 10:50:08.106573933 +0000
1190@@ -242,6 +242,76 @@
1191 large (MPFR_EMAX_MAX);
1192 }
1193
1194+/* Cases where the function overflows on n = 0 when rounding is like
1195+ away from zero. */
1196+static void
1197+overflow0 (mpfr_exp_t emax)
1198+{
1199+ mpfr_exp_t old_emax;
1200+ mpfr_t x, y1, y2;
1201+ int neg, r, op;
1202+ static char *sop[4] = { "mul_2ui", "mul_2si", "div_2ui", "div_2si" };
1203+
1204+ old_emax = mpfr_get_emax ();
1205+ set_emax (emax);
1206+
1207+ mpfr_init2 (x, 8);
1208+ mpfr_inits2 (6, y1, y2, (mpfr_ptr) 0);
1209+
1210+ mpfr_set_inf (x, 1);
1211+ mpfr_nextbelow (x);
1212+
1213+ for (neg = 0; neg <= 1; neg++)
1214+ {
1215+ RND_LOOP (r)
1216+ {
1217+ int inex1, inex2;
1218+ unsigned int flags1, flags2;
1219+
1220+ /* Even if there isn't an overflow (rounding ~ toward zero),
1221+ the result is the same as the one of an overflow. */
1222+ inex1 = mpfr_overflow (y1, (mpfr_rnd_t) r, neg ? -1 : 1);
1223+ flags1 = MPFR_FLAGS_INEXACT;
1224+ if (mpfr_inf_p (y1))
1225+ flags1 |= MPFR_FLAGS_OVERFLOW;
1226+ for (op = 0; op < 4; op++)
1227+ {
1228+ mpfr_clear_flags ();
1229+ inex2 =
1230+ op == 0 ? mpfr_mul_2ui (y2, x, 0, (mpfr_rnd_t) r) :
1231+ op == 1 ? mpfr_mul_2si (y2, x, 0, (mpfr_rnd_t) r) :
1232+ op == 2 ? mpfr_div_2ui (y2, x, 0, (mpfr_rnd_t) r) :
1233+ op == 3 ? mpfr_div_2si (y2, x, 0, (mpfr_rnd_t) r) :
1234+ (MPFR_ASSERTN (0), 0);
1235+ flags2 = __gmpfr_flags;
1236+ if (!(mpfr_equal_p (y1, y2) &&
1237+ SAME_SIGN (inex1, inex2) &&
1238+ flags1 == flags2))
1239+ {
1240+ printf ("Error in overflow0 for %s, mpfr_%s, emax = %"
1241+ MPFR_EXP_FSPEC "d,\nx = ",
1242+ mpfr_print_rnd_mode ((mpfr_rnd_t) r), sop[op],
1243+ (mpfr_eexp_t) emax);
1244+ mpfr_dump (x);
1245+ printf ("Expected ");
1246+ mpfr_dump (y1);
1247+ printf (" with inex = %d, flags =", inex1);
1248+ flags_out (flags1);
1249+ printf ("Got ");
1250+ mpfr_dump (y2);
1251+ printf (" with inex = %d, flags =", inex2);
1252+ flags_out (flags2);
1253+ exit (1);
1254+ }
1255+ }
1256+ }
1257+ mpfr_neg (x, x, MPFR_RNDN);
1258+ }
1259+
1260+ mpfr_clears (x, y1, y2, (mpfr_ptr) 0);
1261+ set_emax (old_emax);
1262+}
1263+
1264 int
1265 main (int argc, char *argv[])
1266 {
1267@@ -334,6 +404,11 @@
1268 underflow0 ();
1269 large0 ();
1270
1271+ if (mpfr_get_emax () != MPFR_EMAX_MAX)
1272+ overflow0 (mpfr_get_emax ());
1273+ overflow0 (MPFR_EMAX_MAX);
1274+ overflow0 (-1);
1275+
1276 tests_end_mpfr ();
1277 return 0;
1278 }
1279diff -Naurd mpfr-3.1.3-a/PATCHES mpfr-3.1.3-b/PATCHES
1280--- mpfr-3.1.3-a/PATCHES 2015-07-17 08:54:48.592799981 +0000
1281+++ mpfr-3.1.3-b/PATCHES 2015-07-17 08:54:48.616811495 +0000
1282@@ -0,0 +1 @@
1283+muldiv-2exp-underflow
1284diff -Naurd mpfr-3.1.3-a/VERSION mpfr-3.1.3-b/VERSION
1285--- mpfr-3.1.3-a/VERSION 2015-07-02 10:50:08.126574142 +0000
1286+++ mpfr-3.1.3-b/VERSION 2015-07-17 08:54:48.616811495 +0000
1287@@ -1 +1 @@
1288-3.1.3-p2
1289+3.1.3-p3
1290diff -Naurd mpfr-3.1.3-a/src/div_2si.c mpfr-3.1.3-b/src/div_2si.c
1291--- mpfr-3.1.3-a/src/div_2si.c 2015-07-02 10:50:08.106573933 +0000
1292+++ mpfr-3.1.3-b/src/div_2si.c 2015-07-17 08:54:48.608807656 +0000
1293@@ -45,7 +45,8 @@
1294 if (rnd_mode == MPFR_RNDN &&
1295 (__gmpfr_emin > MPFR_EMAX_MAX - (n - 1) ||
1296 exp < __gmpfr_emin + (n - 1) ||
1297- (inexact >= 0 && mpfr_powerof2_raw (y))))
1298+ ((MPFR_IS_NEG (y) ? inexact <= 0 : inexact >= 0) &&
1299+ mpfr_powerof2_raw (y))))
1300 rnd_mode = MPFR_RNDZ;
1301 return mpfr_underflow (y, rnd_mode, MPFR_SIGN(y));
1302 }
1303diff -Naurd mpfr-3.1.3-a/src/div_2ui.c mpfr-3.1.3-b/src/div_2ui.c
1304--- mpfr-3.1.3-a/src/div_2ui.c 2015-07-02 10:50:08.106573933 +0000
1305+++ mpfr-3.1.3-b/src/div_2ui.c 2015-07-17 08:54:48.608807656 +0000
1306@@ -44,7 +44,9 @@
1307 if (MPFR_UNLIKELY (n >= diffexp)) /* exp - n <= emin - 1 */
1308 {
1309 if (rnd_mode == MPFR_RNDN &&
1310- (n > diffexp || (inexact >= 0 && mpfr_powerof2_raw (y))))
1311+ (n > diffexp ||
1312+ ((MPFR_IS_NEG (y) ? inexact <= 0 : inexact >= 0) &&
1313+ mpfr_powerof2_raw (y))))
1314 rnd_mode = MPFR_RNDZ;
1315 return mpfr_underflow (y, rnd_mode, MPFR_SIGN (y));
1316 }
1317diff -Naurd mpfr-3.1.3-a/src/mpfr.h mpfr-3.1.3-b/src/mpfr.h
1318--- mpfr-3.1.3-a/src/mpfr.h 2015-07-02 10:50:08.126574142 +0000
1319+++ mpfr-3.1.3-b/src/mpfr.h 2015-07-17 08:54:48.616811495 +0000
1320@@ -27,7 +27,7 @@
1321 #define MPFR_VERSION_MAJOR 3
1322 #define MPFR_VERSION_MINOR 1
1323 #define MPFR_VERSION_PATCHLEVEL 3
1324-#define MPFR_VERSION_STRING "3.1.3-p2"
1325+#define MPFR_VERSION_STRING "3.1.3-p3"
1326
1327 /* Macros dealing with MPFR VERSION */
1328 #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c))
1329diff -Naurd mpfr-3.1.3-a/src/mul_2si.c mpfr-3.1.3-b/src/mul_2si.c
1330--- mpfr-3.1.3-a/src/mul_2si.c 2015-07-02 10:50:08.106573933 +0000
1331+++ mpfr-3.1.3-b/src/mul_2si.c 2015-07-17 08:54:48.608807656 +0000
1332@@ -48,7 +48,8 @@
1333 if (rnd_mode == MPFR_RNDN &&
1334 (__gmpfr_emin > MPFR_EMAX_MAX + (n + 1) ||
1335 exp < __gmpfr_emin - (n + 1) ||
1336- (inexact >= 0 && mpfr_powerof2_raw (y))))
1337+ ((MPFR_IS_NEG (y) ? inexact <= 0 : inexact >= 0) &&
1338+ mpfr_powerof2_raw (y))))
1339 rnd_mode = MPFR_RNDZ;
1340 return mpfr_underflow (y, rnd_mode, MPFR_SIGN(y));
1341 }
1342diff -Naurd mpfr-3.1.3-a/src/version.c mpfr-3.1.3-b/src/version.c
1343--- mpfr-3.1.3-a/src/version.c 2015-07-02 10:50:08.126574142 +0000
1344+++ mpfr-3.1.3-b/src/version.c 2015-07-17 08:54:48.616811495 +0000
1345@@ -25,5 +25,5 @@
1346 const char *
1347 mpfr_get_version (void)
1348 {
1349- return "3.1.3-p2";
1350+ return "3.1.3-p3";
1351 }
1352diff -Naurd mpfr-3.1.3-a/tests/tmul_2exp.c mpfr-3.1.3-b/tests/tmul_2exp.c
1353--- mpfr-3.1.3-a/tests/tmul_2exp.c 2015-07-02 10:50:08.106573933 +0000
1354+++ mpfr-3.1.3-b/tests/tmul_2exp.c 2015-07-17 08:54:48.608807656 +0000
1355@@ -50,77 +50,82 @@
1356 {
1357 mpfr_t x, y, z1, z2;
1358 mpfr_exp_t emin;
1359- int i, k;
1360+ int i, k, s;
1361 int prec;
1362 int rnd;
1363 int div;
1364 int inex1, inex2;
1365 unsigned int flags1, flags2;
1366
1367- /* Test mul_2si(x, e - k), div_2si(x, k - e) and div_2ui(x, k - e)
1368- * with emin = e, x = 1 + i/16, i in { -1, 0, 1 }, and k = 1 to 4,
1369- * by comparing the result with the one of a simple division.
1370+ /* Test mul_2si(x, e - k), div_2si(x, k - e) and div_2ui(x, k - e) with
1371+ * emin = e, x = s * (1 + i/16), i in { -1, 0, 1 }, s in { -1, 1 }, and
1372+ * k = 1 to 4, by comparing the result with the one of a simple division.
1373 */
1374 emin = mpfr_get_emin ();
1375 set_emin (e);
1376 mpfr_inits2 (8, x, y, (mpfr_ptr) 0);
1377 for (i = 15; i <= 17; i++)
1378- {
1379- inex1 = mpfr_set_ui_2exp (x, i, -4, MPFR_RNDN);
1380- MPFR_ASSERTN (inex1 == 0);
1381- for (prec = 6; prec >= 3; prec -= 3)
1382- {
1383- mpfr_inits2 (prec, z1, z2, (mpfr_ptr) 0);
1384- RND_LOOP (rnd)
1385- for (k = 1; k <= 4; k++)
1386- {
1387- /* The following one is assumed to be correct. */
1388- inex1 = mpfr_mul_2si (y, x, e, MPFR_RNDN);
1389- MPFR_ASSERTN (inex1 == 0);
1390- inex1 = mpfr_set_ui (z1, 1 << k, MPFR_RNDN);
1391- MPFR_ASSERTN (inex1 == 0);
1392- mpfr_clear_flags ();
1393- /* Do not use mpfr_div_ui to avoid the optimization
1394- by mpfr_div_2si. */
1395- inex1 = mpfr_div (z1, y, z1, (mpfr_rnd_t) rnd);
1396- flags1 = __gmpfr_flags;
1397-
1398- for (div = 0; div <= 2; div++)
1399+ for (s = 1; s >= -1; s -= 2)
1400+ {
1401+ inex1 = mpfr_set_si_2exp (x, s * i, -4, MPFR_RNDN);
1402+ MPFR_ASSERTN (inex1 == 0);
1403+ for (prec = 6; prec >= 3; prec -= 3)
1404+ {
1405+ mpfr_inits2 (prec, z1, z2, (mpfr_ptr) 0);
1406+ RND_LOOP (rnd)
1407+ for (k = 1; k <= 4; k++)
1408 {
1409+ /* The following one is assumed to be correct. */
1410+ inex1 = mpfr_mul_2si (y, x, e, MPFR_RNDN);
1411+ MPFR_ASSERTN (inex1 == 0);
1412+ inex1 = mpfr_set_ui (z1, 1 << k, MPFR_RNDN);
1413+ MPFR_ASSERTN (inex1 == 0);
1414 mpfr_clear_flags ();
1415- inex2 = div == 0 ?
1416- mpfr_mul_2si (z2, x, e - k, (mpfr_rnd_t) rnd) : div == 1 ?
1417- mpfr_div_2si (z2, x, k - e, (mpfr_rnd_t) rnd) :
1418- mpfr_div_2ui (z2, x, k - e, (mpfr_rnd_t) rnd);
1419- flags2 = __gmpfr_flags;
1420- if (flags1 == flags2 && SAME_SIGN (inex1, inex2) &&
1421- mpfr_equal_p (z1, z2))
1422- continue;
1423- printf ("Error in underflow(");
1424- if (e == MPFR_EMIN_MIN)
1425- printf ("MPFR_EMIN_MIN");
1426- else if (e == emin)
1427- printf ("default emin");
1428- else if (e >= LONG_MIN)
1429- printf ("%ld", (long) e);
1430- else
1431- printf ("<LONG_MIN");
1432- printf (") with mpfr_%s,\nx = %d/16, prec = %d, k = %d, "
1433- "%s\n", div == 0 ? "mul_2si" : div == 1 ?
1434- "div_2si" : "div_2ui", i, prec, k,
1435- mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
1436- printf ("Expected ");
1437- mpfr_out_str (stdout, 16, 0, z1, MPFR_RNDN);
1438- printf (", inex = %d, flags = %u\n", SIGN (inex1), flags1);
1439- printf ("Got ");
1440- mpfr_out_str (stdout, 16, 0, z2, MPFR_RNDN);
1441- printf (", inex = %d, flags = %u\n", SIGN (inex2), flags2);
1442- exit (1);
1443- } /* div */
1444- } /* k */
1445- mpfr_clears (z1, z2, (mpfr_ptr) 0);
1446- } /* prec */
1447- } /* i */
1448+ /* Do not use mpfr_div_ui to avoid the optimization
1449+ by mpfr_div_2si. */
1450+ inex1 = mpfr_div (z1, y, z1, (mpfr_rnd_t) rnd);
1451+ flags1 = __gmpfr_flags;
1452+
1453+ for (div = 0; div <= 2; div++)
1454+ {
1455+ mpfr_clear_flags ();
1456+ inex2 =
1457+ div == 0 ?
1458+ mpfr_mul_2si (z2, x, e - k, (mpfr_rnd_t) rnd) :
1459+ div == 1 ?
1460+ mpfr_div_2si (z2, x, k - e, (mpfr_rnd_t) rnd) :
1461+ mpfr_div_2ui (z2, x, k - e, (mpfr_rnd_t) rnd);
1462+ flags2 = __gmpfr_flags;
1463+ if (flags1 == flags2 && SAME_SIGN (inex1, inex2) &&
1464+ mpfr_equal_p (z1, z2))
1465+ continue;
1466+ printf ("Error in underflow(");
1467+ if (e == MPFR_EMIN_MIN)
1468+ printf ("MPFR_EMIN_MIN");
1469+ else if (e == emin)
1470+ printf ("default emin");
1471+ else if (e >= LONG_MIN)
1472+ printf ("%ld", (long) e);
1473+ else
1474+ printf ("<LONG_MIN");
1475+ printf (") with mpfr_%s,\nx = %d/16, prec = %d, k = %d,"
1476+ " %s\n", div == 0 ? "mul_2si" : div == 1 ?
1477+ "div_2si" : "div_2ui", s * i, prec, k,
1478+ mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
1479+ printf ("Expected ");
1480+ mpfr_out_str (stdout, 16, 0, z1, MPFR_RNDN);
1481+ printf (", inex = %d, flags = %u\n",
1482+ SIGN (inex1), flags1);
1483+ printf ("Got ");
1484+ mpfr_out_str (stdout, 16, 0, z2, MPFR_RNDN);
1485+ printf (", inex = %d, flags = %u\n",
1486+ SIGN (inex2), flags2);
1487+ exit (1);
1488+ } /* div */
1489+ } /* k */
1490+ mpfr_clears (z1, z2, (mpfr_ptr) 0);
1491+ } /* prec */
1492+ } /* i */
1493 mpfr_clears (x, y, (mpfr_ptr) 0);
1494 set_emin (emin);
1495 }
1496diff -Naurd mpfr-3.1.3-a/PATCHES mpfr-3.1.3-b/PATCHES
1497--- mpfr-3.1.3-a/PATCHES 2015-07-17 08:58:21.094987384 +0000
1498+++ mpfr-3.1.3-b/PATCHES 2015-07-17 08:58:21.118986898 +0000
1499@@ -0,0 +1 @@
1500+frexp
1501diff -Naurd mpfr-3.1.3-a/VERSION mpfr-3.1.3-b/VERSION
1502--- mpfr-3.1.3-a/VERSION 2015-07-17 08:54:48.616811495 +0000
1503+++ mpfr-3.1.3-b/VERSION 2015-07-17 08:58:21.118986898 +0000
1504@@ -1 +1 @@
1505-3.1.3-p3
1506+3.1.3-p4
1507diff -Naurd mpfr-3.1.3-a/src/frexp.c mpfr-3.1.3-b/src/frexp.c
1508--- mpfr-3.1.3-a/src/frexp.c 2015-06-19 19:55:09.000000000 +0000
1509+++ mpfr-3.1.3-b/src/frexp.c 2015-07-17 08:58:21.106987142 +0000
1510@@ -26,6 +26,13 @@
1511 mpfr_frexp (mpfr_exp_t *exp, mpfr_ptr y, mpfr_srcptr x, mpfr_rnd_t rnd)
1512 {
1513 int inex;
1514+ unsigned int saved_flags = __gmpfr_flags;
1515+ MPFR_BLOCK_DECL (flags);
1516+
1517+ MPFR_LOG_FUNC
1518+ (("x[%Pu]=%.*Rg rnd=%d", mpfr_get_prec (x), mpfr_log_prec, x, rnd),
1519+ ("y[%Pu]=%.*Rg exp=%" MPFR_EXP_FSPEC "d inex=%d", mpfr_get_prec (y),
1520+ mpfr_log_prec, y, (mpfr_eexp_t) *exp, inex));
1521
1522 if (MPFR_UNLIKELY(MPFR_IS_SINGULAR(x)))
1523 {
1524@@ -49,8 +56,32 @@
1525 }
1526 }
1527
1528- inex = mpfr_set (y, x, rnd);
1529+ MPFR_BLOCK (flags, inex = mpfr_set (y, x, rnd));
1530+ __gmpfr_flags = saved_flags;
1531+
1532+ /* Possible overflow due to the rounding, no possible underflow. */
1533+
1534+ if (MPFR_UNLIKELY (MPFR_OVERFLOW (flags)))
1535+ {
1536+ int inex2;
1537+
1538+ /* An overflow here means that the exponent of y would be larger than
1539+ the one of x, thus x would be rounded to the next power of 2, and
1540+ the returned y should be 1/2 in absolute value, rounded (i.e. with
1541+ possible underflow or overflow). This also implies that x and y are
1542+ different objects, so that the exponent of x has not been lost. */
1543+ MPFR_LOG_MSG (("Internal overflow\n", 0));
1544+ MPFR_ASSERTD (x != y);
1545+ *exp = MPFR_GET_EXP (x) + 1;
1546+ inex2 = mpfr_set_si_2exp (y, MPFR_INT_SIGN (x), -1, rnd);
1547+ MPFR_LOG_MSG (("inex=%d inex2=%d\n", inex, inex2));
1548+ if (inex2 != 0)
1549+ inex = inex2;
1550+ MPFR_RET (inex);
1551+ }
1552+
1553 *exp = MPFR_GET_EXP (y);
1554- MPFR_SET_EXP (y, 0);
1555+ /* Do not use MPFR_SET_EXP because the range has not been checked yet. */
1556+ MPFR_EXP (y) = 0;
1557 return mpfr_check_range (y, inex, rnd);
1558 }
1559diff -Naurd mpfr-3.1.3-a/src/mpfr.h mpfr-3.1.3-b/src/mpfr.h
1560--- mpfr-3.1.3-a/src/mpfr.h 2015-07-17 08:54:48.616811495 +0000
1561+++ mpfr-3.1.3-b/src/mpfr.h 2015-07-17 08:58:21.114986979 +0000
1562@@ -27,7 +27,7 @@
1563 #define MPFR_VERSION_MAJOR 3
1564 #define MPFR_VERSION_MINOR 1
1565 #define MPFR_VERSION_PATCHLEVEL 3
1566-#define MPFR_VERSION_STRING "3.1.3-p3"
1567+#define MPFR_VERSION_STRING "3.1.3-p4"
1568
1569 /* Macros dealing with MPFR VERSION */
1570 #define MPFR_VERSION_NUM(a,b,c) (((a) << 16L) | ((b) << 8) | (c))
1571diff -Naurd mpfr-3.1.3-a/src/version.c mpfr-3.1.3-b/src/version.c
1572--- mpfr-3.1.3-a/src/version.c 2015-07-17 08:54:48.616811495 +0000
1573+++ mpfr-3.1.3-b/src/version.c 2015-07-17 08:58:21.118986898 +0000
1574@@ -25,5 +25,5 @@
1575 const char *
1576 mpfr_get_version (void)
1577 {
1578- return "3.1.3-p3";
1579+ return "3.1.3-p4";
1580 }
1581diff -Naurd mpfr-3.1.3-a/tests/tfrexp.c mpfr-3.1.3-b/tests/tfrexp.c
1582--- mpfr-3.1.3-a/tests/tfrexp.c 2015-06-19 19:55:10.000000000 +0000
1583+++ mpfr-3.1.3-b/tests/tfrexp.c 2015-07-17 08:58:21.106987142 +0000
1584@@ -129,12 +129,115 @@
1585 mpfr_clear (x);
1586 }
1587
1588+static void check1 (void)
1589+{
1590+ mpfr_exp_t emin, emax, e;
1591+ mpfr_t x, y1, y2;
1592+ int r, neg, red;
1593+
1594+ emin = mpfr_get_emin ();
1595+ emax = mpfr_get_emax ();
1596+ set_emin (MPFR_EMIN_MIN);
1597+ set_emax (MPFR_EMAX_MAX);
1598+
1599+ mpfr_init2 (x, 7);
1600+ mpfr_inits2 (4, y1, y2, (mpfr_ptr) 0);
1601+
1602+ mpfr_set_ui_2exp (x, 1, -2, MPFR_RNDN);
1603+ while (mpfr_regular_p (x))
1604+ {
1605+ /* Test the exponents up to 3 and with the maximum exponent
1606+ (to check potential intermediate overflow). */
1607+ if (MPFR_GET_EXP (x) == 4)
1608+ mpfr_set_exp (x, MPFR_EMAX_MAX);
1609+ e = MPFR_GET_EXP (x);
1610+ for (neg = 0; neg < 2; neg++)
1611+ {
1612+ RND_LOOP (r)
1613+ {
1614+ int inex1, inex2;
1615+ mpfr_exp_t e1, e2;
1616+ unsigned int flags1, flags2;
1617+
1618+ for (red = 0; red < 2; red++)
1619+ {
1620+ if (red)
1621+ {
1622+ /* e1: exponent of the rounded value of x. */
1623+ MPFR_ASSERTN (e1 == e || e1 == e + 1);
1624+ set_emin (e);
1625+ set_emax (e);
1626+ mpfr_clear_flags ();
1627+ inex1 = e1 < 0 ?
1628+ mpfr_mul_2ui (y1, x, -e1, (mpfr_rnd_t) r) :
1629+ mpfr_div_2ui (y1, x, e1, (mpfr_rnd_t) r);
1630+ flags1 = __gmpfr_flags;
1631+ }
1632+ else
1633+ {
1634+ inex1 = mpfr_set (y1, x, (mpfr_rnd_t) r);
1635+ e1 = MPFR_IS_INF (y1) ? e + 1 : MPFR_GET_EXP (y1);
1636+ flags1 = inex1 != 0 ? MPFR_FLAGS_INEXACT : 0;
1637+ }
1638+ mpfr_clear_flags ();
1639+ inex2 = mpfr_frexp (&e2, y2, x, (mpfr_rnd_t) r);
1640+ flags2 = __gmpfr_flags;
1641+ set_emin (MPFR_EMIN_MIN);
1642+ set_emax (MPFR_EMAX_MAX);
1643+ if ((!red || e == 0) &&
1644+ (! mpfr_regular_p (y2) || MPFR_GET_EXP (y2) != 0))
1645+ {
1646+ printf ("Error in check1 for %s, red = %d, x = ",
1647+ mpfr_print_rnd_mode ((mpfr_rnd_t) r), red);
1648+ mpfr_dump (x);
1649+ printf ("Expected 1/2 <= |y| < 1, got y = ");
1650+ mpfr_dump (y2);
1651+ exit (1);
1652+ }
1653+ if (!red)
1654+ {
1655+ if (e2 > 0)
1656+ mpfr_mul_2ui (y2, y2, e2, MPFR_RNDN);
1657+ else if (e2 < 0)
1658+ mpfr_div_2ui (y2, y2, -e2, MPFR_RNDN);
1659+ }
1660+ if (! (SAME_SIGN (inex1, inex2) &&
1661+ mpfr_equal_p (y1, y2) &&
1662+ flags1 == flags2))
1663+ {
1664+ printf ("Error in check1 for %s, red = %d, x = ",
1665+ mpfr_print_rnd_mode ((mpfr_rnd_t) r), red);
1666+ mpfr_dump (x);
1667+ printf ("Expected y1 = ");
1668+ mpfr_dump (y1);
1669+ printf ("Got y2 = ");
1670+ mpfr_dump (y2);
1671+ printf ("Expected inex ~= %d, got %d\n", inex1, inex2);
1672+ printf ("Expected flags:");
1673+ flags_out (flags1);
1674+ printf ("Got flags: ");
1675+ flags_out (flags2);
1676+ exit (1);
1677+ }
1678+ }
1679+ }
1680+ mpfr_neg (x, x, MPFR_RNDN);
1681+ }
1682+ mpfr_nextabove (x);
1683+ }
1684+
1685+ mpfr_clears (x, y1, y2, (mpfr_ptr) 0);
1686+ set_emin (emin);
1687+ set_emax (emax);
1688+}
1689+
1690 int
1691 main (int argc, char *argv[])
1692 {
1693 tests_start_mpfr ();
1694
1695 check_special ();
1696+ check1 ();
1697
1698 tests_end_mpfr ();
1699 return 0;