00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #ifndef __MMX_INTERVAL_HPP
00014 #define __MMX_INTERVAL_HPP
00015 #include <basix/function.hpp>
00016 #include <numerix/rounded.hpp>
00017 namespace mmx {
00018 #define TMPL_DEF template<typename C,typename V=std_interval>
00019 #define TMPL template<typename C,typename V>
00020 #define Interval interval<C,V>
00021 #define Up V::template rnd_helper<C>::UV
00022 #define Down V::template rnd_helper<C>::DV
00023 TMPL class interval;
00024 TMPL inline C lower (const Interval& z);
00025 TMPL inline C upper (const Interval& z);
00026
00027
00028
00029
00030
00031 struct std_interval {
00032 template<typename C>
00033 struct rnd_helper {
00034 typedef Round_up(C) UV;
00035 typedef Round_down(C) DV;
00036 };
00037 };
00038
00039
00040
00041
00042
00043 TMPL_DEF
00044 class interval {
00045 MMX_ALLOCATORS
00046 private:
00047 C l;
00048 C r;
00049
00050 public:
00051 typedef C boundary_type;
00052
00053 inline interval (): l (0), r (0) {}
00054 template<typename T> inline interval (const T& x):
00055 l (Down::template as<C> (x)),
00056 r (Up::template as<C> (x)) {}
00057 template<typename T> inline interval (const interval<T,V>& x):
00058 l (Down::template as<C> (lower (x))),
00059 r (Up::template as<C> (upper (x))) {}
00060 template<typename LT,typename UT> inline
00061 interval (const LT& l2, const UT& u2):
00062 l (Down::template as<C> (l2)),
00063 r (Up::template as<C> (u2)) {}
00064 inline interval (const C& x): l (x), r (x) {}
00065 inline interval (const C& l2, const C& r2): l (l2), r (r2) {}
00066 friend C lower LESSGTR (const Interval& x);
00067 friend C upper LESSGTR (const Interval& x);
00068
00069 static inline Interval nan ();
00070 static inline Interval maximal ();
00071 static inline Interval minimal ();
00072 static inline Interval infinity ();
00073 static inline Interval fuzzy ();
00074 static inline Interval log2 ();
00075 static inline Interval pi ();
00076 static inline Interval euler ();
00077
00078 inline Interval& operator += (const Interval& x);
00079 inline Interval& operator -= (const Interval& x);
00080 inline Interval& operator <<= (const xint& shift);
00081 inline Interval& operator >>= (const xint& shift);
00082 };
00083 DEFINE_UNARY_FORMAT_1 (interval)
00084
00085 UNARY_RETURN_TYPE(TMPL,center_op,Interval,C);
00086 UNARY_RETURN_TYPE(TMPL,radius_op,Interval,C);
00087 TMPL struct rounding_helper<Interval >: public rounding_helper<C> {};
00088
00089
00090
00091
00092
00093 template<typename C, typename V, typename C2, typename R2>
00094 struct make_ball_helper<Interval,C2,R2> {
00095 static inline Interval val (const C2& c, const R2& r) {
00096 typedef Round_up(C2) Rnd_Up;
00097 typedef Round_down(C2) Rnd_Down;
00098 return Interval (Rnd_Down::template sub_as<C2> (c, r),
00099 Rnd_Up::template add_as<C2> (c, r));
00100 }
00101 };
00102
00103 template<typename C, typename V, typename C2>
00104 struct make_interval_helper<Interval,C2> {
00105 static inline Interval val (const C2& l, const C2& r) {
00106 return Interval (l, r); }
00107 };
00108
00109 TMPL inline C lower (const Interval& x) { return x.l; }
00110 TMPL inline C upper (const Interval& x) { return x.r; }
00111 TMPL inline C center (const Interval& x) {
00112 return decexp2 (lower (x) + upper (x)); }
00113 TMPL inline C radius (const Interval& x) {
00114 if (is_infinite (x)) return 0;
00115 return Up::add (decexp2 (Up::sub (upper (x), lower (x))),
00116 multiplicative_error (center (x))); }
00117 TMPL inline format<C> CF (const Interval& x) {
00118 return get_format (lower (x)); }
00119
00120 TMPL inline Interval abs (const Interval& x);
00121 TMPL inline C abs_down (const Interval& x) { return lower (abs (x)); }
00122 TMPL inline C abs_up (const Interval& x) { return upper (abs (x)); }
00123 TMPL inline C bnd_down (const Interval& x) { return lower (x); }
00124 TMPL inline C bnd_up (const Interval& x) { return upper (x); }
00125
00126
00127
00128
00129
00130 TMPL inline Interval copy (const Interval& x) {
00131 return Interval (copy (lower (x)), copy (upper (x))); }
00132 TMPL inline Interval duplicate (const Interval& x) {
00133 return Interval (duplicate (lower (x)), duplicate (upper (x))); }
00134
00135 TMPL inline bool is_exact_zero (const Interval& x) {
00136 return is_exact_zero (lower (x)) && is_exact_zero (upper (x)); }
00137
00138 template<typename Op, typename C, typename V> nat
00139 unary_hash (const Interval& x) {
00140 nat h= Op::op (lower (x));
00141 return (h<<3) ^ (h<<15) ^ (h>>29) ^ Op::op (upper (x));
00142 }
00143
00144 template<typename Op, typename C, typename V> bool
00145 binary_test (const Interval& x1, const Interval& x2) {
00146 return Op::op (lower (x1), lower (x2)) && Op::op (upper (x1), upper (x2));
00147 }
00148
00149 template<typename S1, typename D> interval<D>
00150 map (const function_1<D,Argument(S1) >& fun, const interval<S1>& x) {
00151 return interval<D> (fun (lower (x)), fun (upper (x)));
00152 }
00153
00154 template<typename S1, typename D> interval<D>
00155 map (D (*fun) (const S1&), const interval<S1>& x) {
00156 return interval<D> (fun (lower (x)), fun (upper (x)));
00157 }
00158
00159 EXACT_IDENTITY_OP_SUGAR(TMPL,Interval)
00160 HARD_IDENTITY_OP_SUGAR(TMPL,Interval)
00161 TMPL nat hash (const Interval& x) { return 12721; }
00162
00163 TMPL syntactic
00164 flatten (const Interval& x) {
00165 if (is_nan (x)) return Nan (syntactic);
00166 if (lower (x) == 0 && upper (x) == 0) return 0;
00167 return flatten_range (lower (x), upper (x));
00168
00169 }
00170
00171 TMPL syntactic
00172 flatten (const Interval& x, xnat dd) {
00173 if (is_nan (x)) return Nan (syntactic);
00174 if (lower (x) == 0 && upper (x) == 0) return 0;
00175 return flatten_range (lower (x), upper (x), dd);
00176 }
00177
00178 TMPL
00179 struct binary_helper<Interval >: public void_binary_helper<Interval > {
00180 static inline string short_type_name () {
00181 return "I" * Short_type_name (C); }
00182 static inline generic full_type_name () {
00183 return gen ("Interval", Full_type_name (C)); }
00184 static inline generic disassemble (const Interval& x) {
00185 return gen_vec (as<generic> (lower (x)), as<generic> (upper (x))); }
00186 static inline Interval assemble (const generic& v) {
00187 return Interval (as<C> (vector_access (v, 0)),
00188 as<C> (vector_access (v, 1))); }
00189 static inline void write (const port& out, const Interval& x) {
00190 binary_write<C> (out, lower (x));
00191 binary_write<C> (out, upper (x)); }
00192 static inline Interval read (const port& in) {
00193 C l= binary_read<C> (in);
00194 C r= binary_read<C> (in);
00195 return Interval (l, r); }
00196 };
00197
00198
00199
00200
00201
00202 TMPL inline bool
00203 intersect (Interval& r, const Interval &a, const Interval &b) {
00204 if ((upper (a) < lower (b)) || (lower (a) > upper (b))) return false;
00205 r = Interval (max (lower (a), lower (b)), min (upper (a), upper (b)));
00206 return true;
00207 }
00208
00209 TMPL inline bool
00210 is_zero (const Interval& x) {
00211 return is_finite (x) && sign (lower (x)) * sign (upper (x)) <= 0;
00212 }
00213
00214 TMPL inline bool
00215 is_non_zero (const Interval& x) {
00216 return !is_finite (x) || sign (lower (x)) * sign (upper (x)) > 0;
00217 }
00218
00219 TMPL inline bool
00220 operator == (const Interval& x1, const Interval& x2) {
00221 if (is_nan (x1) || is_nan (x2)) return is_nan (x1) && is_nan (x2);
00222 return lower (x1) <= upper (x2) && lower (x2) <= upper (x1);
00223 }
00224
00225 TMPL inline bool
00226 operator != (const Interval& x1, const Interval& x2) {
00227 if (is_nan (x1) || is_nan (x2)) return !(is_nan (x1) && is_nan (x2));
00228 return lower (x1) > upper (x2) || lower (x2) > upper (x1);
00229 }
00230
00231 TMPL inline bool
00232 operator <= (const Interval& x1, const Interval& x2) {
00233 if (is_nan (x1) || is_nan (x2)) return false;
00234 return lower (x1) <= upper (x2);
00235 }
00236
00237 TMPL inline bool
00238 operator < (const Interval& x1, const Interval& x2) {
00239 if (is_nan (x1) || is_nan (x2)) return false;
00240 return lower (x2) > upper (x1);
00241 }
00242
00243 TMPL inline bool
00244 operator >= (const Interval& x1, const Interval& x2) {
00245 if (is_nan (x1) || is_nan (x2)) return false;
00246 return lower (x2) <= upper (x1);
00247 }
00248
00249 TMPL inline bool
00250 operator > (const Interval& x1, const Interval& x2) {
00251 if (is_nan (x1) || is_nan (x2)) return false;
00252 return lower (x1) > upper (x2);
00253 }
00254
00255 EQUAL_INT_SUGAR(TMPL,Interval);
00256 COMPARE_INT_SUGAR(TMPL,Interval);
00257 EQUAL_SCALAR_SUGAR(TMPL,Interval,C);
00258 COMPARE_SCALAR_SUGAR(TMPL,Interval,C);
00259 EQUAL_SCALAR_SUGAR_BIS(TMPL,Interval,C);
00260 COMPARE_SCALAR_SUGAR_BIS(TMPL,Interval,C);
00261
00262
00263
00264
00265
00266 TMPL inline Interval
00267 abs (const Interval& x) {
00268 if (is_nan (x)) return x;
00269 C r= max (abs (lower (x)), abs (upper (x)));
00270 if (sign (lower (x)) * sign (upper (x)) <= 0) return Interval (C(0), r);
00271 return Interval (min (abs (lower (x)), abs (upper (x))), r);
00272 }
00273
00274 TMPL Interval
00275 min (const Interval& x1, const Interval& x2) {
00276 return Interval (min (lower (x1), lower (x2)), min (upper (x1), upper (x2)));
00277 }
00278
00279 TMPL Interval
00280 max (const Interval& x1, const Interval& x2) {
00281 return Interval (max (lower (x1), lower (x2)), max (upper (x1), upper (x2)));
00282 }
00283
00284 TMPL Interval
00285 floor (const Interval& x) {
00286 return Interval (Down::floor (lower (x)), Up::floor (upper (x)));
00287 }
00288
00289 TMPL Interval
00290 trunc (const Interval& x) {
00291 return Interval (Down::trunc (lower (x)), Up::trunc (upper (x)));
00292 }
00293
00294 TMPL Interval
00295 ceil (const Interval& x) {
00296 return Interval (Down::ceil (lower (x)), Up::ceil (upper (x)));
00297 }
00298
00299 TMPL Interval
00300 round (const Interval& x) {
00301 return Interval (Down::round (lower (x)), Up::round (upper (x)));
00302 }
00303
00304
00305
00306
00307
00308 TMPL inline void set_default (Interval& x) {
00309 x= Interval (Default (C), Default (C)); }
00310 TMPL inline void set_nan (Interval& x) {
00311 x= Interval (Nan (C), Nan (C)); }
00312 TMPL inline void set_maximal (Interval& x) {
00313 x= Interval (Maximal (C), Maximal (C)); }
00314 TMPL inline void set_minimal (Interval& x) {
00315 x= Interval (Minimal (C), Minimal (C)); }
00316 TMPL inline void set_infinity (Interval& x) {
00317 x= Interval (Infinity (C), Infinity (C)); }
00318 TMPL inline void set_fuzz (Interval& x) {
00319 x= Interval (Minimal (C), Maximal (C)); }
00320 TMPL inline void set_smallest (Interval& x) {
00321 x= Interval (Smallest (C), Smallest (C)); }
00322 TMPL inline void set_largest (Interval& x) {
00323 x= Interval (Largest (C), Largest (C)); }
00324 TMPL inline void set_accuracy (Interval& x) {
00325 x= Interval (Accuracy (C), Accuracy (C)); }
00326 TMPL inline void set_log2 (Interval& x) {
00327 x= Interval (Down::template log2<C>(), Up::template log2<C>()); }
00328 TMPL inline void set_pi (Interval& x) {
00329 x= Interval (Down::template pi<C>(), Up::template pi<C>()); }
00330 TMPL inline void set_euler (Interval& x) {
00331 x= Interval (Down::template euler<C>(), Up::template euler<C>()); }
00332 TMPL inline void set_catalan (Interval& x) {
00333 x= Interval (Down::template catalan<C>(), Up::template catalan<C>()); }
00334 TMPL inline Interval times_infinity (const Interval& x) {
00335 return x * infinity_cst<Interval > (); }
00336
00337
00338
00339
00340
00341 TMPL inline Interval
00342 operator - (const Interval& x) {
00343 return Interval (-upper (x), -lower (x));
00344 }
00345
00346 TMPL inline Interval
00347 operator + (const Interval& x1, const Interval& x2) {
00348 return Interval (Down::add (lower (x1), lower (x2)),
00349 Up::add (upper (x1), upper (x2)));
00350 }
00351
00352 TMPL inline Interval
00353 operator + (const Interval& x1, const C& x2) {
00354 return Interval (Down::add (lower (x1), x2), Up::add (upper (x1), x2));
00355 }
00356
00357 TMPL inline Interval
00358 operator + (const C& x1, const Interval& x2) {
00359 return Interval (Down::add (x1, lower (x2)), Up::add (x1, upper (x2)));
00360 }
00361
00362 TMPL inline Interval
00363 operator - (const Interval& x1, const Interval& x2) {
00364 return Interval (Down::sub (lower (x1), upper (x2)),
00365 Up::sub (upper (x1), lower (x2)));
00366 }
00367
00368 TMPL inline Interval
00369 operator - (const Interval& x1, const C& x2) {
00370 return Interval (Down::add (lower (x1), x2), Up::add (upper (x1), x2));
00371 }
00372
00373 TMPL inline Interval
00374 operator - (const C& x1, const Interval& x2) {
00375 return Interval (Down::add (x1, upper (x2)), Up::add (x1, lower (x2)));
00376 }
00377
00378 TMPL Interval
00379 operator * (const Interval& x1, const Interval& x2) {
00380 if (lower (x1) >= 0)
00381 return Interval (lower (x2) >= 0?
00382 Down::mul (lower (x1), lower (x2)):
00383 Down::mul (upper (x1), lower (x2)),
00384 upper (x2) >= 0?
00385 Up::mul (upper (x1), upper (x2)):
00386 Up::mul (lower (x1), upper (x2)));
00387 else if (upper (x1) <= 0)
00388 return Interval (upper (x2) >= 0?
00389 Down::mul (lower (x1), upper (x2)):
00390 Down::mul (upper (x1), upper (x2)),
00391 lower (x2) >= 0?
00392 Up::mul (upper (x1), lower (x2)):
00393 Up::mul (lower (x1), lower (x2)));
00394 else if (lower (x2) >= 0)
00395 return Interval (Down::mul (lower (x1), upper (x2)),
00396 Up::mul (upper (x1), upper (x2)));
00397 else if (upper (x2) <= 0)
00398 return Interval (Down::mul (upper (x1), lower (x2)),
00399 Up::mul (lower (x1), lower (x2)));
00400 else
00401 return Interval (min (Down::mul (lower (x1), upper (x2)),
00402 Down::mul (upper (x1), lower (x2))),
00403 max (Up::mul (lower (x1), lower (x2)),
00404 Up::mul (upper (x1), upper (x2))));
00405 }
00406
00407 TMPL inline Interval operator * (const Interval& x1, const C& x2) {
00408 return x1 * Interval (x2); }
00409 TMPL inline Interval operator * (const C& x1, const Interval& x2) {
00410 return Interval (x1) * x2; }
00411
00412 TMPL inline Interval
00413 square (const Interval& x) {
00414 if (lower (x) >= 0)
00415 return Interval (Down::mul (lower (x), lower (x)),
00416 Up::mul (upper (x), upper (x)));
00417 else if (upper (x) <= 0)
00418 return Interval (Down::mul (upper (x), upper (x)),
00419 Up::mul (lower (x), lower (x)));
00420 else return Interval (Down::mul (lower (x), upper (x)),
00421 max (Up::mul (lower (x), lower (x)),
00422 Up::mul (upper (x), upper (x))));
00423 }
00424
00425 TMPL Interval
00426 operator / (const Interval& x1, const Interval& x2) {
00427 if (is_zero (x2)) return Nan (Interval);
00428 else if (lower (x2) > 0)
00429 return Interval (lower (x1) > 0?
00430 Down::div (lower (x1), upper (x2)):
00431 Down::div (lower (x1), lower (x2)),
00432 upper (x1) > 0?
00433 Up::div (upper (x1), lower (x2)):
00434 Up::div (upper (x1), upper (x2)));
00435 else if (upper (x2) < 0)
00436 return Interval (upper (x1) > 0?
00437 Down::div (upper (x1), upper (x2)):
00438 Down::div (upper (x1), lower (x2)),
00439 lower (x1) > 0?
00440 Up::div (lower (x1), lower (x2)):
00441 Up::div (lower (x1), upper (x2)));
00442 else return Nan (Interval);
00443 }
00444
00445 TMPL inline Interval operator / (const Interval& x1, const C& x2) {
00446 return x1 / Interval (x2); }
00447 TMPL inline Interval operator / (const C& x1, const Interval& x2) {
00448 return Interval (x1) / x2; }
00449
00450 ARITH_INT_SUGAR(TMPL,Interval);
00451
00452
00453
00454
00455
00456 TMPL inline Interval&
00457 Interval::operator += (const Interval& x) {
00458 l= Down::add (l, lower (x));
00459 r= Up::add (r, upper (x));
00460 return *this;
00461 }
00462
00463 TMPL inline Interval&
00464 Interval::operator -= (const Interval& x) {
00465 l= Down::sub (l, upper (x));
00466 r= Up::sub (r, lower (x));
00467 return *this;
00468 }
00469
00470
00471
00472
00473
00474 TMPL Interval
00475 sqrt (const Interval& x) {
00476 return Interval (Down::sqrt (lower (x)), Up::sqrt (upper (x)));
00477 }
00478
00479 TMPL Interval
00480 hypot (const Interval& x, const Interval& y) {
00481 Interval z= square (x) + square (y);
00482 if (z > 0) return sqrt (z);
00483 else return Interval (0, Up::sqrt (upper (z)));
00484 }
00485
00486 TMPL Interval
00487 exp (const Interval& x) {
00488 return Interval (Down::exp (lower (x)), Up::exp (upper (x)));
00489 }
00490
00491 TMPL Interval
00492 log (const Interval& x) {
00493 return Interval (Down::log (lower (x)), Up::log (upper (x)));
00494 }
00495
00496 TMPL inline Interval
00497 pow (const Interval& x, const Interval& y) {
00498 return exp (y * log (x));
00499 }
00500
00501 TMPL Interval
00502 cos (const Interval& x) {
00503 C pi_lo= promote (314159, lower (x)) / promote (100000, lower (x));
00504 if (Up::sub (upper (x), lower (x)) >= pi_lo)
00505 return Interval (promote (-1, lower (x)), promote (1, lower (x)));
00506 else {
00507 int s1= -sign (sin (lower (x)));
00508 int s2= -sign (sin (upper (x)));
00509 if (s1 >= 0 && s2 >= 0)
00510 return Interval (Down::cos (lower (x)), Up::cos (upper (x)));
00511 else if (s1 <= 0 && s2 <= 0)
00512 return Interval (Down::cos (upper (x)), Down::cos (lower (x)));
00513 else if (s1 >= 0)
00514 return Interval (min (Down::cos (lower (x)), Down::cos (upper (x))),
00515 promote (1, lower (x)));
00516 else
00517 return Interval (promote (-1, lower (x)),
00518 max (Down::cos (lower (x)), Down::cos (upper (x))));
00519 }
00520 }
00521
00522 TMPL Interval
00523 sin (const Interval& x) {
00524 C pi_lo= promote (314159, lower (x)) / promote (100000, lower (x));
00525 if (Up::sub (upper (x), lower (x)) >= pi_lo)
00526 return Interval (promote (-1, lower (x)), promote (1, lower (x)));
00527 else {
00528 int s1= sign (cos (lower (x)));
00529 int s2= sign (cos (upper (x)));
00530 if (s1 >= 0 && s2 >= 0)
00531 return Interval (Down::sin (lower (x)), Up::sin (upper (x)));
00532 else if (s1 <= 0 && s2 <= 0)
00533 return Interval (Down::sin (upper (x)), Down::sin (lower (x)));
00534 else if (s1 >= 0)
00535 return Interval (min (Down::sin (lower (x)), Down::sin (upper (x))),
00536 promote (1, lower (x)));
00537 else
00538 return Interval (promote (-1, lower (x)),
00539 max (Down::sin (lower (x)), Down::sin (upper (x))));
00540 }
00541 }
00542
00543 TMPL Interval
00544 tan (const Interval& x) {
00545 return sin (x) / cos (x);
00546 }
00547
00548 TMPL inline Interval
00549 cosh (const Interval& x) {
00550 if (is_nan (x)) return Nan (Interval);
00551 else if (lower (x) >= 0)
00552 return Interval (Down::cosh (lower (x)), Up::cosh (upper (x)));
00553 else if (upper (x) <= 0)
00554 return Interval (Down::cosh (upper (x)), Up::cosh (lower (x)));
00555 else return Interval (1, max (Up::cosh (lower (x)), Up::cosh (upper (x))));
00556 }
00557
00558 TMPL inline Interval
00559 sinh (const Interval& x) {
00560 return Interval (Down::sinh (lower (x)), Up::sinh (upper (x)));
00561 }
00562
00563 TMPL Interval
00564 tanh (const Interval& x) {
00565 return sinh (x) / cosh (x);
00566 }
00567
00568 TMPL Interval
00569 acos (const Interval& x) {
00570 return Interval (Down::acos (upper (x)), Up::acos (lower (x)));
00571 }
00572
00573 TMPL Interval
00574 asin (const Interval& x) {
00575 return Interval (Down::asin (lower (x)), Up::asin (upper (x)));
00576 }
00577
00578 TMPL Interval
00579 atan (const Interval& x) {
00580 return Interval (Down::atan (lower (x)), Up::atan (upper (x)));
00581 }
00582
00583 TMPL Interval
00584 atan2 (const Interval& y, const Interval& x) {
00585 if (sign (lower (x)) > 0) return atan (y/x);
00586 if (sign (lower (y)) > 0) return atan (-x/y) + (Pi (Interval) >> 1);
00587 if (sign (upper (y)) < 0) return atan (-x/y) - (Pi (Interval) >> 1);
00588 if (sign (upper (x)) < 0) return atan (y/x) + Pi (Interval);
00589 C pi= Up::template pi<C> ();
00590 return Interval (-pi, pi);
00591 }
00592
00593 INV_TRIGO_SUGAR(TMPL,Interval)
00594 INV_HYPER_SUGAR(TMPL,Interval)
00595 ARG_HYPER_SUGAR(TMPL,Interval)
00596
00597
00598
00599
00600
00601 TMPL inline bool is_finite (const Interval& x) {
00602 return is_finite (lower (x)) && is_finite (upper (x)); }
00603 TMPL inline bool is_infinite (const Interval& x) {
00604 return is_infinite (lower (x)) && upper (x) == lower (x); }
00605 TMPL inline bool is_fuzz (const Interval& x) {
00606 return is_infinite (lower (x)) && is_infinite (upper (x)) &&
00607 lower (x) != upper (x); }
00608 TMPL inline bool is_nan (const Interval& x) {
00609 return is_nan (lower (x)) || is_nan (upper (x)); }
00610 TMPL inline bool is_reliable (const Interval& x) {
00611 (void) x; return true; }
00612
00613 TMPL inline Interval change_precision (const Interval& x, xnat prec) {
00614 return Interval (change_precision (lower (x), prec),
00615 change_precision (upper (x), prec)); }
00616 TMPL inline xnat precision (const Interval& x) {
00617 return min (precision (lower (x)), precision (upper (x))); }
00618
00619 TMPL inline xint exponent (const Interval& x) {
00620 return max (exponent (lower (x)), exponent (upper (x))); }
00621 TMPL inline double magnitude (const Interval& x) {
00622 return max (magnitude (lower (x)), magnitude (upper (x))); }
00623 TMPL inline Interval operator << (const Interval& x, const xint& shift) {
00624 return Interval (incexp2 (lower (x), shift), incexp2 (upper (x), shift)); }
00625 TMPL inline Interval operator >> (const Interval& x, const xint& shift) {
00626 return Interval (decexp2 (lower (x), shift), decexp2 (upper (x), shift)); }
00627
00628 TMPL inline Interval&
00629 Interval::operator <<= (const xint& shift) {
00630 incexp2_assign (l, shift);
00631 incexp2_assign (r, shift);
00632 return *this;
00633 }
00634
00635 TMPL inline Interval&
00636 Interval::operator >>= (const xint& shift) {
00637 decexp2_assign (l, shift);
00638 decexp2_assign (r, shift);
00639 return *this;
00640 }
00641
00642 TMPL inline Interval sharpen (const Interval& x) {
00643 return Interval (center (x)); }
00644 TMPL inline Interval blur (const Interval& x, const C& r) {
00645 return x + Interval (-r, r); }
00646 TMPL inline Interval blur (const Interval& x, const Interval& r) {
00647 return blur (x, abs_up (r)); }
00648
00649 #undef TMPL_DEF
00650 #undef TMPL
00651 #undef Interval
00652 #undef Up
00653 #undef Down
00654 }
00655 #endif // __MMX_INTERVAL_HPP