1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/os/ossrv/ossrv_pub/boost_apis/boost/numeric/interval/arith.hpp Fri Jun 15 03:10:57 2012 +0200
1.3 @@ -0,0 +1,305 @@
1.4 +/* Boost interval/arith.hpp template implementation file
1.5 + *
1.6 + * Copyright 2000 Jens Maurer
1.7 + * Copyright 2002-2003 Hervé Brönnimann, Guillaume Melquiond, Sylvain Pion
1.8 + *
1.9 + * Distributed under the Boost Software License, Version 1.0.
1.10 + * (See accompanying file LICENSE_1_0.txt or
1.11 + * copy at http://www.boost.org/LICENSE_1_0.txt)
1.12 + */
1.13 +
1.14 +#ifndef BOOST_NUMERIC_INTERVAL_ARITH_HPP
1.15 +#define BOOST_NUMERIC_INTERVAL_ARITH_HPP
1.16 +
1.17 +#include <boost/config.hpp>
1.18 +#include <boost/numeric/interval/interval.hpp>
1.19 +#include <boost/numeric/interval/detail/bugs.hpp>
1.20 +#include <boost/numeric/interval/detail/test_input.hpp>
1.21 +#include <boost/numeric/interval/detail/division.hpp>
1.22 +#include <algorithm>
1.23 +
1.24 +namespace boost {
1.25 +namespace numeric {
1.26 +
1.27 +/*
1.28 + * Basic arithmetic operators
1.29 + */
1.30 +
1.31 +template<class T, class Policies> inline
1.32 +const interval<T, Policies>& operator+(const interval<T, Policies>& x)
1.33 +{
1.34 + return x;
1.35 +}
1.36 +
1.37 +template<class T, class Policies> inline
1.38 +interval<T, Policies> operator-(const interval<T, Policies>& x)
1.39 +{
1.40 + if (interval_lib::detail::test_input(x))
1.41 + return interval<T, Policies>::empty();
1.42 + return interval<T, Policies>(-x.upper(), -x.lower(), true);
1.43 +}
1.44 +
1.45 +template<class T, class Policies> inline
1.46 +interval<T, Policies>& interval<T, Policies>::operator+=(const interval<T, Policies>& r)
1.47 +{
1.48 + if (interval_lib::detail::test_input(*this, r))
1.49 + set_empty();
1.50 + else {
1.51 + typename Policies::rounding rnd;
1.52 + set(rnd.add_down(low, r.low), rnd.add_up(up, r.up));
1.53 + }
1.54 + return *this;
1.55 +}
1.56 +
1.57 +template<class T, class Policies> inline
1.58 +interval<T, Policies>& interval<T, Policies>::operator+=(const T& r)
1.59 +{
1.60 + if (interval_lib::detail::test_input(*this, r))
1.61 + set_empty();
1.62 + else {
1.63 + typename Policies::rounding rnd;
1.64 + set(rnd.add_down(low, r), rnd.add_up(up, r));
1.65 + }
1.66 + return *this;
1.67 +}
1.68 +
1.69 +template<class T, class Policies> inline
1.70 +interval<T, Policies>& interval<T, Policies>::operator-=(const interval<T, Policies>& r)
1.71 +{
1.72 + if (interval_lib::detail::test_input(*this, r))
1.73 + set_empty();
1.74 + else {
1.75 + typename Policies::rounding rnd;
1.76 + set(rnd.sub_down(low, r.up), rnd.sub_up(up, r.low));
1.77 + }
1.78 + return *this;
1.79 +}
1.80 +
1.81 +template<class T, class Policies> inline
1.82 +interval<T, Policies>& interval<T, Policies>::operator-=(const T& r)
1.83 +{
1.84 + if (interval_lib::detail::test_input(*this, r))
1.85 + set_empty();
1.86 + else {
1.87 + typename Policies::rounding rnd;
1.88 + set(rnd.sub_down(low, r), rnd.sub_up(up, r));
1.89 + }
1.90 + return *this;
1.91 +}
1.92 +
1.93 +template<class T, class Policies> inline
1.94 +interval<T, Policies>& interval<T, Policies>::operator*=(const interval<T, Policies>& r)
1.95 +{
1.96 + return *this = *this * r;
1.97 +}
1.98 +
1.99 +template<class T, class Policies> inline
1.100 +interval<T, Policies>& interval<T, Policies>::operator*=(const T& r)
1.101 +{
1.102 + return *this = r * *this;
1.103 +}
1.104 +
1.105 +template<class T, class Policies> inline
1.106 +interval<T, Policies>& interval<T, Policies>::operator/=(const interval<T, Policies>& r)
1.107 +{
1.108 + return *this = *this / r;
1.109 +}
1.110 +
1.111 +template<class T, class Policies> inline
1.112 +interval<T, Policies>& interval<T, Policies>::operator/=(const T& r)
1.113 +{
1.114 + return *this = *this / r;
1.115 +}
1.116 +
1.117 +template<class T, class Policies> inline
1.118 +interval<T, Policies> operator+(const interval<T, Policies>& x,
1.119 + const interval<T, Policies>& y)
1.120 +{
1.121 + if (interval_lib::detail::test_input(x, y))
1.122 + return interval<T, Policies>::empty();
1.123 + typename Policies::rounding rnd;
1.124 + return interval<T,Policies>(rnd.add_down(x.lower(), y.lower()),
1.125 + rnd.add_up (x.upper(), y.upper()), true);
1.126 +}
1.127 +
1.128 +template<class T, class Policies> inline
1.129 +interval<T, Policies> operator+(const T& x, const interval<T, Policies>& y)
1.130 +{
1.131 + if (interval_lib::detail::test_input(x, y))
1.132 + return interval<T, Policies>::empty();
1.133 + typename Policies::rounding rnd;
1.134 + return interval<T,Policies>(rnd.add_down(x, y.lower()),
1.135 + rnd.add_up (x, y.upper()), true);
1.136 +}
1.137 +
1.138 +template<class T, class Policies> inline
1.139 +interval<T, Policies> operator+(const interval<T, Policies>& x, const T& y)
1.140 +{ return y + x; }
1.141 +
1.142 +template<class T, class Policies> inline
1.143 +interval<T, Policies> operator-(const interval<T, Policies>& x,
1.144 + const interval<T, Policies>& y)
1.145 +{
1.146 + if (interval_lib::detail::test_input(x, y))
1.147 + return interval<T, Policies>::empty();
1.148 + typename Policies::rounding rnd;
1.149 + return interval<T,Policies>(rnd.sub_down(x.lower(), y.upper()),
1.150 + rnd.sub_up (x.upper(), y.lower()), true);
1.151 +}
1.152 +
1.153 +template<class T, class Policies> inline
1.154 +interval<T, Policies> operator-(const T& x, const interval<T, Policies>& y)
1.155 +{
1.156 + if (interval_lib::detail::test_input(x, y))
1.157 + return interval<T, Policies>::empty();
1.158 + typename Policies::rounding rnd;
1.159 + return interval<T,Policies>(rnd.sub_down(x, y.upper()),
1.160 + rnd.sub_up (x, y.lower()), true);
1.161 +}
1.162 +
1.163 +template<class T, class Policies> inline
1.164 +interval<T, Policies> operator-(const interval<T, Policies>& x, const T& y)
1.165 +{
1.166 + if (interval_lib::detail::test_input(x, y))
1.167 + return interval<T, Policies>::empty();
1.168 + typename Policies::rounding rnd;
1.169 + return interval<T,Policies>(rnd.sub_down(x.lower(), y),
1.170 + rnd.sub_up (x.upper(), y), true);
1.171 +}
1.172 +
1.173 +template<class T, class Policies> inline
1.174 +interval<T, Policies> operator*(const interval<T, Policies>& x,
1.175 + const interval<T, Policies>& y)
1.176 +{
1.177 + BOOST_USING_STD_MIN();
1.178 + BOOST_USING_STD_MAX();
1.179 + typedef interval<T, Policies> I;
1.180 + if (interval_lib::detail::test_input(x, y))
1.181 + return I::empty();
1.182 + typename Policies::rounding rnd;
1.183 + const T& xl = x.lower();
1.184 + const T& xu = x.upper();
1.185 + const T& yl = y.lower();
1.186 + const T& yu = y.upper();
1.187 +
1.188 + if (interval_lib::user::is_neg(xl))
1.189 + if (interval_lib::user::is_pos(xu))
1.190 + if (interval_lib::user::is_neg(yl))
1.191 + if (interval_lib::user::is_pos(yu)) // M * M
1.192 + return I(min BOOST_PREVENT_MACRO_SUBSTITUTION(rnd.mul_down(xl, yu), rnd.mul_down(xu, yl)),
1.193 + max BOOST_PREVENT_MACRO_SUBSTITUTION(rnd.mul_up (xl, yl), rnd.mul_up (xu, yu)), true);
1.194 + else // M * N
1.195 + return I(rnd.mul_down(xu, yl), rnd.mul_up(xl, yl), true);
1.196 + else
1.197 + if (interval_lib::user::is_pos(yu)) // M * P
1.198 + return I(rnd.mul_down(xl, yu), rnd.mul_up(xu, yu), true);
1.199 + else // M * Z
1.200 + return I(static_cast<T>(0), static_cast<T>(0), true);
1.201 + else
1.202 + if (interval_lib::user::is_neg(yl))
1.203 + if (interval_lib::user::is_pos(yu)) // N * M
1.204 + return I(rnd.mul_down(xl, yu), rnd.mul_up(xl, yl), true);
1.205 + else // N * N
1.206 + return I(rnd.mul_down(xu, yu), rnd.mul_up(xl, yl), true);
1.207 + else
1.208 + if (interval_lib::user::is_pos(yu)) // N * P
1.209 + return I(rnd.mul_down(xl, yu), rnd.mul_up(xu, yl), true);
1.210 + else // N * Z
1.211 + return I(static_cast<T>(0), static_cast<T>(0), true);
1.212 + else
1.213 + if (interval_lib::user::is_pos(xu))
1.214 + if (interval_lib::user::is_neg(yl))
1.215 + if (interval_lib::user::is_pos(yu)) // P * M
1.216 + return I(rnd.mul_down(xu, yl), rnd.mul_up(xu, yu), true);
1.217 + else // P * N
1.218 + return I(rnd.mul_down(xu, yl), rnd.mul_up(xl, yu), true);
1.219 + else
1.220 + if (interval_lib::user::is_pos(yu)) // P * P
1.221 + return I(rnd.mul_down(xl, yl), rnd.mul_up(xu, yu), true);
1.222 + else // P * Z
1.223 + return I(static_cast<T>(0), static_cast<T>(0), true);
1.224 + else // Z * ?
1.225 + return I(static_cast<T>(0), static_cast<T>(0), true);
1.226 +}
1.227 +
1.228 +template<class T, class Policies> inline
1.229 +interval<T, Policies> operator*(const T& x, const interval<T, Policies>& y)
1.230 +{
1.231 + typedef interval<T, Policies> I;
1.232 + if (interval_lib::detail::test_input(x, y))
1.233 + return I::empty();
1.234 + typename Policies::rounding rnd;
1.235 + const T& yl = y.lower();
1.236 + const T& yu = y.upper();
1.237 + // x is supposed not to be infinite
1.238 + if (interval_lib::user::is_neg(x))
1.239 + return I(rnd.mul_down(x, yu), rnd.mul_up(x, yl), true);
1.240 + else if (interval_lib::user::is_zero(x))
1.241 + return I(static_cast<T>(0), static_cast<T>(0), true);
1.242 + else
1.243 + return I(rnd.mul_down(x, yl), rnd.mul_up(x, yu), true);
1.244 +}
1.245 +
1.246 +template<class T, class Policies> inline
1.247 +interval<T, Policies> operator*(const interval<T, Policies>& x, const T& y)
1.248 +{ return y * x; }
1.249 +
1.250 +template<class T, class Policies> inline
1.251 +interval<T, Policies> operator/(const interval<T, Policies>& x,
1.252 + const interval<T, Policies>& y)
1.253 +{
1.254 + if (interval_lib::detail::test_input(x, y))
1.255 + return interval<T, Policies>::empty();
1.256 + if (zero_in(y))
1.257 + if (!interval_lib::user::is_zero(y.lower()))
1.258 + if (!interval_lib::user::is_zero(y.upper()))
1.259 + return interval_lib::detail::div_zero(x);
1.260 + else
1.261 + return interval_lib::detail::div_negative(x, y.lower());
1.262 + else
1.263 + if (!interval_lib::user::is_zero(y.upper()))
1.264 + return interval_lib::detail::div_positive(x, y.upper());
1.265 + else
1.266 + return interval<T, Policies>::empty();
1.267 + else
1.268 + return interval_lib::detail::div_non_zero(x, y);
1.269 +}
1.270 +
1.271 +template<class T, class Policies> inline
1.272 +interval<T, Policies> operator/(const T& x, const interval<T, Policies>& y)
1.273 +{
1.274 + if (interval_lib::detail::test_input(x, y))
1.275 + return interval<T, Policies>::empty();
1.276 + if (zero_in(y))
1.277 + if (!interval_lib::user::is_zero(y.lower()))
1.278 + if (!interval_lib::user::is_zero(y.upper()))
1.279 + return interval_lib::detail::div_zero<T, Policies>(x);
1.280 + else
1.281 + return interval_lib::detail::div_negative<T, Policies>(x, y.lower());
1.282 + else
1.283 + if (!interval_lib::user::is_zero(y.upper()))
1.284 + return interval_lib::detail::div_positive<T, Policies>(x, y.upper());
1.285 + else
1.286 + return interval<T, Policies>::empty();
1.287 + else
1.288 + return interval_lib::detail::div_non_zero(x, y);
1.289 +}
1.290 +
1.291 +template<class T, class Policies> inline
1.292 +interval<T, Policies> operator/(const interval<T, Policies>& x, const T& y)
1.293 +{
1.294 + if (interval_lib::detail::test_input(x, y) || interval_lib::user::is_zero(y))
1.295 + return interval<T, Policies>::empty();
1.296 + typename Policies::rounding rnd;
1.297 + const T& xl = x.lower();
1.298 + const T& xu = x.upper();
1.299 + if (interval_lib::user::is_neg(y))
1.300 + return interval<T, Policies>(rnd.div_down(xu, y), rnd.div_up(xl, y), true);
1.301 + else
1.302 + return interval<T, Policies>(rnd.div_down(xl, y), rnd.div_up(xu, y), true);
1.303 +}
1.304 +
1.305 +} // namespace numeric
1.306 +} // namespace boost
1.307 +
1.308 +#endif // BOOST_NUMERIC_INTERVAL_ARITH_HPP