1 ///////////////////////////////////////////////////////////////////////////////
2 // Copyright 2011 John Maddock. Distributed under the Boost
3 // Software License, Version 1.0. (See accompanying file
4 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
5
6 #ifndef BOOST_MP_GENERIC_INTERCONVERT_HPP
7 #define BOOST_MP_GENERIC_INTERCONVERT_HPP
8
9 #include <boost/multiprecision/detail/default_ops.hpp>
10
11 #ifdef BOOST_MSVC
12 #pragma warning(push)
13 #pragma warning(disable:4127 6326)
14 #endif
15
16 namespace boost{ namespace multiprecision{ namespace detail{
17
18 template <class To, class From>
do_cast(const From & from)19 inline To do_cast(const From & from)
20 {
21 return static_cast<To>(from);
22 }
23 template <class To, class B, ::boost::multiprecision::expression_template_option et>
do_cast(const number<B,et> & from)24 inline To do_cast(const number<B, et>& from)
25 {
26 return from.template convert_to<To>();
27 }
28
29 template <class To, class From>
generic_interconvert(To & to,const From & from,const mpl::int_<number_kind_floating_point> &,const mpl::int_<number_kind_integer> &)30 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_floating_point>& /*to_type*/, const mpl::int_<number_kind_integer>& /*from_type*/)
31 {
32 using default_ops::eval_get_sign;
33 using default_ops::eval_bitwise_and;
34 using default_ops::eval_convert_to;
35 using default_ops::eval_right_shift;
36 using default_ops::eval_ldexp;
37 using default_ops::eval_add;
38 // smallest unsigned type handled natively by "From" is likely to be it's limb_type:
39 typedef typename canonical<unsigned char, From>::type limb_type;
40 // get the corresponding type that we can assign to "To":
41 typedef typename canonical<limb_type, To>::type to_type;
42 From t(from);
43 bool is_neg = eval_get_sign(t) < 0;
44 if(is_neg)
45 t.negate();
46 // Pick off the first limb:
47 limb_type limb;
48 limb_type mask = ~static_cast<limb_type>(0);
49 From fl;
50 eval_bitwise_and(fl, t, mask);
51 eval_convert_to(&limb, fl);
52 to = static_cast<to_type>(limb);
53 eval_right_shift(t, std::numeric_limits<limb_type>::digits);
54 //
55 // Then keep picking off more limbs until "t" is zero:
56 //
57 To l;
58 unsigned shift = std::numeric_limits<limb_type>::digits;
59 while(!eval_is_zero(t))
60 {
61 eval_bitwise_and(fl, t, mask);
62 eval_convert_to(&limb, fl);
63 l = static_cast<to_type>(limb);
64 eval_right_shift(t, std::numeric_limits<limb_type>::digits);
65 eval_ldexp(l, l, shift);
66 eval_add(to, l);
67 shift += std::numeric_limits<limb_type>::digits;
68 }
69 //
70 // Finish off by setting the sign:
71 //
72 if(is_neg)
73 to.negate();
74 }
75
76 template <class To, class From>
generic_interconvert(To & to,const From & from,const mpl::int_<number_kind_integer> &,const mpl::int_<number_kind_integer> &)77 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_integer>& /*to_type*/, const mpl::int_<number_kind_integer>& /*from_type*/)
78 {
79 using default_ops::eval_get_sign;
80 using default_ops::eval_bitwise_and;
81 using default_ops::eval_convert_to;
82 using default_ops::eval_right_shift;
83 using default_ops::eval_left_shift;
84 using default_ops::eval_bitwise_or;
85 using default_ops::eval_is_zero;
86 // smallest unsigned type handled natively by "From" is likely to be it's limb_type:
87 typedef typename canonical<unsigned char, From>::type limb_type;
88 // get the corresponding type that we can assign to "To":
89 typedef typename canonical<limb_type, To>::type to_type;
90 From t(from);
91 bool is_neg = eval_get_sign(t) < 0;
92 if(is_neg)
93 t.negate();
94 // Pick off the first limb:
95 limb_type limb;
96 limb_type mask = static_cast<limb_type>(~static_cast<limb_type>(0));
97 From fl;
98 eval_bitwise_and(fl, t, mask);
99 eval_convert_to(&limb, fl);
100 to = static_cast<to_type>(limb);
101 eval_right_shift(t, std::numeric_limits<limb_type>::digits);
102 //
103 // Then keep picking off more limbs until "t" is zero:
104 //
105 To l;
106 unsigned shift = std::numeric_limits<limb_type>::digits;
107 while(!eval_is_zero(t))
108 {
109 eval_bitwise_and(fl, t, mask);
110 eval_convert_to(&limb, fl);
111 l = static_cast<to_type>(limb);
112 eval_right_shift(t, std::numeric_limits<limb_type>::digits);
113 eval_left_shift(l, shift);
114 eval_bitwise_or(to, l);
115 shift += std::numeric_limits<limb_type>::digits;
116 }
117 //
118 // Finish off by setting the sign:
119 //
120 if(is_neg)
121 to.negate();
122 }
123
124 template <class To, class From>
generic_interconvert(To & to,const From & from,const mpl::int_<number_kind_floating_point> &,const mpl::int_<number_kind_floating_point> &)125 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_floating_point>& /*to_type*/, const mpl::int_<number_kind_floating_point>& /*from_type*/)
126 {
127 #ifdef BOOST_MSVC
128 #pragma warning(push)
129 #pragma warning(disable:4127)
130 #endif
131 //
132 // The code here only works when the radix of "From" is 2, we could try shifting by other
133 // radixes but it would complicate things.... use a string conversion when the radix is other
134 // than 2:
135 //
136 if(std::numeric_limits<number<From> >::radix != 2)
137 {
138 to = from.str(0, std::ios_base::fmtflags()).c_str();
139 return;
140 }
141
142
143 typedef typename canonical<unsigned char, To>::type ui_type;
144
145 using default_ops::eval_fpclassify;
146 using default_ops::eval_add;
147 using default_ops::eval_subtract;
148 using default_ops::eval_convert_to;
149
150 //
151 // First classify the input, then handle the special cases:
152 //
153 int c = eval_fpclassify(from);
154
155 if(c == (int)FP_ZERO)
156 {
157 to = ui_type(0);
158 return;
159 }
160 else if(c == (int)FP_NAN)
161 {
162 to = static_cast<const char*>("nan");
163 return;
164 }
165 else if(c == (int)FP_INFINITE)
166 {
167 to = static_cast<const char*>("inf");
168 if(eval_get_sign(from) < 0)
169 to.negate();
170 return;
171 }
172
173 typename From::exponent_type e;
174 From f, term;
175 to = ui_type(0);
176
177 eval_frexp(f, from, &e);
178
179 static const int shift = std::numeric_limits<boost::intmax_t>::digits - 1;
180
181 while(!eval_is_zero(f))
182 {
183 // extract int sized bits from f:
184 eval_ldexp(f, f, shift);
185 eval_floor(term, f);
186 e -= shift;
187 eval_ldexp(to, to, shift);
188 typename boost::multiprecision::detail::canonical<boost::intmax_t, To>::type ll;
189 eval_convert_to(&ll, term);
190 eval_add(to, ll);
191 eval_subtract(f, term);
192 }
193 typedef typename To::exponent_type to_exponent;
194 if((e > (std::numeric_limits<to_exponent>::max)()) || (e < (std::numeric_limits<to_exponent>::min)()))
195 {
196 to = static_cast<const char*>("inf");
197 if(eval_get_sign(from) < 0)
198 to.negate();
199 return;
200 }
201 eval_ldexp(to, to, static_cast<to_exponent>(e));
202 #ifdef BOOST_MSVC
203 #pragma warning(pop)
204 #endif
205 }
206
207 template <class To, class From>
generic_interconvert(To & to,const From & from,const mpl::int_<number_kind_rational> &,const mpl::int_<number_kind_rational> &)208 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_rational>& /*to_type*/, const mpl::int_<number_kind_rational>& /*from_type*/)
209 {
210 typedef typename component_type<number<To> >::type to_component_type;
211
212 number<From> t(from);
213 to_component_type n(numerator(t)), d(denominator(t));
214 using default_ops::assign_components;
215 assign_components(to, n.backend(), d.backend());
216 }
217
218 template <class To, class From>
generic_interconvert(To & to,const From & from,const mpl::int_<number_kind_rational> &,const mpl::int_<number_kind_integer> &)219 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_rational>& /*to_type*/, const mpl::int_<number_kind_integer>& /*from_type*/)
220 {
221 typedef typename component_type<number<To> >::type to_component_type;
222
223 number<From> t(from);
224 to_component_type n(t), d(1);
225 using default_ops::assign_components;
226 assign_components(to, n.backend(), d.backend());
227 }
228
229 template <class R, class LargeInteger>
230 R safe_convert_to_float(const LargeInteger& i)
231 {
232 using std::ldexp;
233 if(!i)
234 return R(0);
235 if(std::numeric_limits<R>::is_specialized && std::numeric_limits<R>::max_exponent)
236 {
237 LargeInteger val(i);
238 if(val.sign() < 0)
239 val = -val;
240 unsigned mb = msb(val);
241 if(mb >= std::numeric_limits<R>::max_exponent)
242 {
243 int scale_factor = (int)mb + 1 - std::numeric_limits<R>::max_exponent;
244 BOOST_ASSERT(scale_factor >= 1);
245 val >>= scale_factor;
246 R result = val.template convert_to<R>();
247 if(std::numeric_limits<R>::digits == 0 || std::numeric_limits<R>::digits >= std::numeric_limits<R>::max_exponent)
248 {
249 //
250 // Calculate and add on the remainder, only if there are more
251 // digits in the mantissa that the size of the exponent, in
252 // other words if we are dropping digits in the conversion
253 // otherwise:
254 //
255 LargeInteger remainder(i);
256 remainder &= (LargeInteger(1) << scale_factor) - 1;
257 result += ldexp(safe_convert_to_float<R>(remainder), -scale_factor);
258 }
259 return i.sign() < 0 ? static_cast<R>(-result) : result;
260 }
261 }
262 return i.template convert_to<R>();
263 }
264
265 template <class To, class Integer>
266 inline typename disable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
generic_convert_rational_to_float_imp(To & result,const Integer & n,const Integer & d,const mpl::true_ &)267 generic_convert_rational_to_float_imp(To& result, const Integer& n, const Integer& d, const mpl::true_&)
268 {
269 //
270 // If we get here, then there's something about one type or the other
271 // that prevents an exactly rounded result from being calculated
272 // (or at least it's not clear how to implement such a thing).
273 //
274 using default_ops::eval_divide;
275 number<To> fn(safe_convert_to_float<number<To> >(n)), fd(safe_convert_to_float<number<To> >(d));
276 eval_divide(result, fn.backend(), fd.backend());
277 }
278 template <class To, class Integer>
279 inline typename enable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
generic_convert_rational_to_float_imp(To & result,const Integer & n,const Integer & d,const mpl::true_ &)280 generic_convert_rational_to_float_imp(To& result, const Integer& n, const Integer& d, const mpl::true_&)
281 {
282 //
283 // If we get here, then there's something about one type or the other
284 // that prevents an exactly rounded result from being calculated
285 // (or at least it's not clear how to implement such a thing).
286 //
287 To fd(safe_convert_to_float<To>(d));
288 result = safe_convert_to_float<To>(n);
289 result /= fd;
290 }
291
292 template <class To, class Integer>
293 typename enable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
generic_convert_rational_to_float_imp(To & result,Integer & num,Integer & denom,const mpl::false_ &)294 generic_convert_rational_to_float_imp(To& result, Integer& num, Integer& denom, const mpl::false_&)
295 {
296 //
297 // If we get here, then the precision of type To is known, and the integer type is unbounded
298 // so we can use integer division plus manipulation of the remainder to get an exactly
299 // rounded result.
300 //
301 if(num == 0)
302 {
303 result = 0;
304 return;
305 }
306 bool s = false;
307 if(num < 0)
308 {
309 s = true;
310 num = -num;
311 }
312 int denom_bits = msb(denom);
313 int shift = std::numeric_limits<To>::digits + denom_bits - msb(num) + 1;
314 if(shift > 0)
315 num <<= shift;
316 else if(shift < 0)
317 denom <<= boost::multiprecision::detail::unsigned_abs(shift);
318 Integer q, r;
319 divide_qr(num, denom, q, r);
320 int q_bits = msb(q);
321 if(q_bits == std::numeric_limits<To>::digits)
322 {
323 //
324 // Round up if 2 * r > denom:
325 //
326 r <<= 1;
327 int c = r.compare(denom);
328 if(c > 0)
329 ++q;
330 else if((c == 0) && (q & 1u))
331 {
332 ++q;
333 }
334 }
335 else
336 {
337 BOOST_ASSERT(q_bits == 1 + std::numeric_limits<To>::digits);
338 //
339 // We basically already have the rounding info:
340 //
341 if(q & 1u)
342 {
343 if(r || (q & 2u))
344 ++q;
345 }
346 }
347 using std::ldexp;
348 result = do_cast<To>(q);
349 result = ldexp(result, -shift);
350 if(s)
351 result = -result;
352 }
353 template <class To, class Integer>
354 inline typename disable_if_c<is_number<To>::value || is_floating_point<To>::value>::type
generic_convert_rational_to_float_imp(To & result,Integer & num,Integer & denom,const mpl::false_ & tag)355 generic_convert_rational_to_float_imp(To& result, Integer& num, Integer& denom, const mpl::false_& tag)
356 {
357 number<To> t;
358 generic_convert_rational_to_float_imp(t, num, denom, tag);
359 result = t.backend();
360 }
361
362 template <class To, class From>
generic_convert_rational_to_float(To & result,const From & f)363 inline void generic_convert_rational_to_float(To& result, const From& f)
364 {
365 //
366 // Type From is always a Backend to number<>, or an
367 // instance of number<>, but we allow
368 // To to be either a Backend type, or a real number type,
369 // that way we can call this from generic conversions, and
370 // from specific conversions to built in types.
371 //
372 typedef typename mpl::if_c<is_number<From>::value, From, number<From> >::type actual_from_type;
373 typedef typename mpl::if_c<is_number<To>::value || is_floating_point<To>::value, To, number<To> >::type actual_to_type;
374 typedef typename component_type<actual_from_type>::type integer_type;
375 typedef mpl::bool_<!std::numeric_limits<integer_type>::is_specialized
376 || std::numeric_limits<integer_type>::is_bounded
377 || !std::numeric_limits<actual_to_type>::is_specialized
378 || !std::numeric_limits<actual_to_type>::is_bounded
379 || (std::numeric_limits<actual_to_type>::radix != 2)> dispatch_tag;
380
381 integer_type n(numerator(static_cast<actual_from_type>(f))), d(denominator(static_cast<actual_from_type>(f)));
382 generic_convert_rational_to_float_imp(result, n, d, dispatch_tag());
383 }
384
385 template <class To, class From>
generic_interconvert(To & to,const From & from,const mpl::int_<number_kind_floating_point> &,const mpl::int_<number_kind_rational> &)386 inline void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_floating_point>& /*to_type*/, const mpl::int_<number_kind_rational>& /*from_type*/)
387 {
388 generic_convert_rational_to_float(to, from);
389 }
390
391 template <class To, class From>
generic_interconvert_float2rational(To & to,const From & from,const mpl::int_<2> &)392 void generic_interconvert_float2rational(To& to, const From& from, const mpl::int_<2>& /*radix*/)
393 {
394 typedef typename mpl::front<typename To::unsigned_types>::type ui_type;
395 static const int shift = std::numeric_limits<boost::long_long_type>::digits;
396 typename From::exponent_type e;
397 typename component_type<number<To> >::type num, denom;
398 number<From> val(from);
399 val = frexp(val, &e);
400 while(val)
401 {
402 val = ldexp(val, shift);
403 e -= shift;
404 boost::long_long_type ll = boost::math::lltrunc(val);
405 val -= ll;
406 num <<= shift;
407 num += ll;
408 }
409 denom = ui_type(1u);
410 if(e < 0)
411 denom <<= -e;
412 else if(e > 0)
413 num <<= e;
414 assign_components(to, num.backend(), denom.backend());
415 }
416
417 template <class To, class From, int Radix>
generic_interconvert_float2rational(To & to,const From & from,const mpl::int_<Radix> &)418 void generic_interconvert_float2rational(To& to, const From& from, const mpl::int_<Radix>& /*radix*/)
419 {
420 //
421 // This is almost the same as the binary case above, but we have to use
422 // scalbn and ilogb rather than ldexp and frexp, we also only extract
423 // one Radix digit at a time which is terribly inefficient!
424 //
425 typedef typename mpl::front<typename To::unsigned_types>::type ui_type;
426 typename From::exponent_type e;
427 typename component_type<To>::type num, denom;
428 number<From> val(from);
429 e = ilogb(val);
430 val = scalbn(val, -e);
431 while(val)
432 {
433 boost::long_long_type ll = boost::math::lltrunc(val);
434 val -= ll;
435 val = scalbn(val, 1);
436 num *= Radix;
437 num += ll;
438 --e;
439 }
440 ++e;
441 denom = ui_type(Radix);
442 denom = pow(denom, abs(e));
443 if(e > 0)
444 {
445 num *= denom;
446 denom = 1;
447 }
448 assign_components(to, num, denom);
449 }
450
451 template <class To, class From>
generic_interconvert(To & to,const From & from,const mpl::int_<number_kind_rational> &,const mpl::int_<number_kind_floating_point> &)452 void generic_interconvert(To& to, const From& from, const mpl::int_<number_kind_rational>& /*to_type*/, const mpl::int_<number_kind_floating_point>& /*from_type*/)
453 {
454 generic_interconvert_float2rational(to, from, mpl::int_<std::numeric_limits<number<From> >::radix>());
455 }
456
457 }}} // namespaces
458
459 #ifdef BOOST_MSVC
460 #pragma warning(pop)
461 #endif
462
463 #endif // BOOST_MP_GENERIC_INTERCONVERT_HPP
464
465