interval/detail/division.hpp

0001 /* Boost interval/detail/division.hpp file
0002  *
0003  * Copyright 2003 Guillaume Melquiond, Sylvain Pion
0004  *
0005  * Distributed under the Boost Software License, Version 1.0.
0006  * (See accompanying file LICENSE_1_0.txt or
0007  * copy at http://www.boost.org/LICENSE_1_0.txt)
0008  */
0009
0010 #ifndef BOOST_NUMERIC_INTERVAL_DETAIL_DIVISION_HPP
0011 #define BOOST_NUMERIC_INTERVAL_DETAIL_DIVISION_HPP
0012
0013 #include <boost/numeric/interval/detail/interval_prototype.hpp>
0014 #include <boost/numeric/interval/detail/bugs.hpp>
0015 #include <boost/numeric/interval/detail/test_input.hpp>
0016 #include <boost/numeric/interval/rounded_arith.hpp>
0017 #include <algorithm>
0018
0019 namespace boost {
0020 namespace numeric {
0021 namespace interval_lib {
0022 namespace detail {
0023
0024 template<class T, class Policies> inline
0025 interval<T, Policies> div_non_zero(const interval<T, Policies>& x,
0026                                    const interval<T, Policies>& y)
0027 {
0028   // assert(!in_zero(y));
0029   typename Policies::rounding rnd;
0030   typedef interval<T, Policies> I;
0031   const T& xl = x.lower();
0032   const T& xu = x.upper();
0033   const T& yl = y.lower();
0034   const T& yu = y.upper();
0035   if (::boost::numeric::interval_lib::user::is_neg(xu))
0036     if (::boost::numeric::interval_lib::user::is_neg(yu))
0037       return I(rnd.div_down(xu, yl), rnd.div_up(xl, yu), true);
0038     else
0039       return I(rnd.div_down(xl, yl), rnd.div_up(xu, yu), true);
0040   else if (::boost::numeric::interval_lib::user::is_neg(xl))
0041     if (::boost::numeric::interval_lib::user::is_neg(yu))
0042       return I(rnd.div_down(xu, yu), rnd.div_up(xl, yu), true);
0043     else
0044       return I(rnd.div_down(xl, yl), rnd.div_up(xu, yl), true);
0045   else
0046     if (::boost::numeric::interval_lib::user::is_neg(yu))
0047       return I(rnd.div_down(xu, yu), rnd.div_up(xl, yl), true);
0048     else
0049       return I(rnd.div_down(xl, yu), rnd.div_up(xu, yl), true);
0050 }
0051
0052 template<class T, class Policies> inline
0053 interval<T, Policies> div_non_zero(const T& x, const interval<T, Policies>& y)
0054 {
0055   // assert(!in_zero(y));
0056   typename Policies::rounding rnd;
0057   typedef interval<T, Policies> I;
0058   const T& yl = y.lower();
0059   const T& yu = y.upper();
0060   if (::boost::numeric::interval_lib::user::is_neg(x))
0061     return I(rnd.div_down(x, yl), rnd.div_up(x, yu), true);
0062   else
0063     return I(rnd.div_down(x, yu), rnd.div_up(x, yl), true);
0064 }
0065
0066 template<class T, class Policies> inline
0067 interval<T, Policies> div_positive(const interval<T, Policies>& x, const T& yu)
0068 {
0069   // assert(::boost::numeric::interval_lib::user::is_pos(yu));
0070   if (::boost::numeric::interval_lib::user::is_zero(x.lower()) &&
0071       ::boost::numeric::interval_lib::user::is_zero(x.upper()))
0072     return x;
0073   typename Policies::rounding rnd;
0074   typedef interval<T, Policies> I;
0075   const T& xl = x.lower();
0076   const T& xu = x.upper();
0077   typedef typename Policies::checking checking;
0078   if (::boost::numeric::interval_lib::user::is_neg(xu))
0079     return I(checking::neg_inf(), rnd.div_up(xu, yu), true);
0080   else if (::boost::numeric::interval_lib::user::is_neg(xl))
0081     return I(checking::neg_inf(), checking::pos_inf(), true);
0082   else
0083     return I(rnd.div_down(xl, yu), checking::pos_inf(), true);
0084 }
0085
0086 template<class T, class Policies> inline
0087 interval<T, Policies> div_positive(const T& x, const T& yu)
0088 {
0089   // assert(::boost::numeric::interval_lib::user::is_pos(yu));
0090   typedef interval<T, Policies> I;
0091   if (::boost::numeric::interval_lib::user::is_zero(x))
0092     return I(static_cast<T>(0), static_cast<T>(0), true);
0093   typename Policies::rounding rnd;
0094   typedef typename Policies::checking checking;
0095   if (::boost::numeric::interval_lib::user::is_neg(x))
0096     return I(checking::neg_inf(), rnd.div_up(x, yu), true);
0097   else
0098     return I(rnd.div_down(x, yu), checking::pos_inf(), true);
0099 }
0100
0101 template<class T, class Policies> inline
0102 interval<T, Policies> div_negative(const interval<T, Policies>& x, const T& yl)
0103 {
0104   // assert(::boost::numeric::interval_lib::user::is_neg(yl));
0105   if (::boost::numeric::interval_lib::user::is_zero(x.lower()) &&
0106       ::boost::numeric::interval_lib::user::is_zero(x.upper()))
0107     return x;
0108   typename Policies::rounding rnd;
0109   typedef interval<T, Policies> I;
0110   const T& xl = x.lower();
0111   const T& xu = x.upper();
0112   typedef typename Policies::checking checking;
0113   if (::boost::numeric::interval_lib::user::is_neg(xu))
0114     return I(rnd.div_down(xu, yl), checking::pos_inf(), true);
0115   else if (::boost::numeric::interval_lib::user::is_neg(xl))
0116     return I(checking::neg_inf(), checking::pos_inf(), true);
0117   else
0118     return I(checking::neg_inf(), rnd.div_up(xl, yl), true);
0119 }
0120
0121 template<class T, class Policies> inline
0122 interval<T, Policies> div_negative(const T& x, const T& yl)
0123 {
0124   // assert(::boost::numeric::interval_lib::user::is_neg(yl));
0125   typedef interval<T, Policies> I;
0126   if (::boost::numeric::interval_lib::user::is_zero(x))
0127     return I(static_cast<T>(0), static_cast<T>(0), true);
0128   typename Policies::rounding rnd;
0129   typedef typename Policies::checking checking;
0130   if (::boost::numeric::interval_lib::user::is_neg(x))
0131     return I(rnd.div_down(x, yl), checking::pos_inf(), true);
0132   else
0133     return I(checking::neg_inf(), rnd.div_up(x, yl), true);
0134 }
0135
0136 template<class T, class Policies> inline
0137 interval<T, Policies> div_zero(const interval<T, Policies>& x)
0138 {
0139   if (::boost::numeric::interval_lib::user::is_zero(x.lower()) &&
0140       ::boost::numeric::interval_lib::user::is_zero(x.upper()))
0141     return x;
0142   else return interval<T, Policies>::whole();
0143 }
0144
0145 template<class T, class Policies> inline
0146 interval<T, Policies> div_zero(const T& x)
0147 {
0148   if (::boost::numeric::interval_lib::user::is_zero(x))
0149     return interval<T, Policies>(static_cast<T>(0), static_cast<T>(0), true);
0150   else return interval<T, Policies>::whole();
0151 }
0152
0153 template<class T, class Policies> inline
0154 interval<T, Policies> div_zero_part1(const interval<T, Policies>& x,
0155                                      const interval<T, Policies>& y, bool& b)
0156 {
0157   // assert(::boost::numeric::interval_lib::user::is_neg(y.lower()) && ::boost::numeric::interval_lib::user::is_pos(y.upper()));
0158   if (::boost::numeric::interval_lib::user::is_zero(x.lower()) && ::boost::numeric::interval_lib::user::is_zero(x.upper()))
0159     { b = false; return x; }
0160   typename Policies::rounding rnd;
0161   typedef interval<T, Policies> I;
0162   const T& xl = x.lower();
0163   const T& xu = x.upper();
0164   const T& yl = y.lower();
0165   const T& yu = y.upper();
0166   typedef typename Policies::checking checking;
0167   if (::boost::numeric::interval_lib::user::is_neg(xu))
0168     { b = true;  return I(checking::neg_inf(), rnd.div_up(xu, yu), true); }
0169   else if (::boost::numeric::interval_lib::user::is_neg(xl))
0170     { b = false; return I(checking::neg_inf(), checking::pos_inf(), true); }
0171   else
0172     { b = true;  return I(checking::neg_inf(), rnd.div_up(xl, yl), true); }
0173 }
0174
0175 template<class T, class Policies> inline
0176 interval<T, Policies> div_zero_part2(const interval<T, Policies>& x,
0177                                      const interval<T, Policies>& y)
0178 {
0179   // assert(::boost::numeric::interval_lib::user::is_neg(y.lower()) && ::boost::numeric::interval_lib::user::is_pos(y.upper()) && (div_zero_part1(x, y, b), b));
0180   typename Policies::rounding rnd;
0181   typedef interval<T, Policies> I;
0182   typedef typename Policies::checking checking;
0183   if (::boost::numeric::interval_lib::user::is_neg(x.upper()))
0184     return I(rnd.div_down(x.upper(), y.lower()), checking::pos_inf(), true);
0185   else
0186     return I(rnd.div_down(x.lower(), y.upper()), checking::pos_inf(), true);
0187 }
0188
0189 } // namespace detail
0190 } // namespace interval_lib
0191 } // namespace numeric
0192 } // namespace boost
0193
0194 #endif // BOOST_NUMERIC_INTERVAL_DETAIL_DIVISION_HPP