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