COCONUT API: prointerval.h Source File

00001 // Projective Interval Class -*- C++ -*-
00002 
00003 // $Id: prointerval.h 847 2008-12-03 23:36:05Z hermannb $
00004 // Copyright (C) 2001-2004 Hermann Schichl and Xuan-Ha VU
00005 //
00006 // This file is part of the COCONUT API.  This library
00007 // is free software; you can redistribute it and/or modify it under the
00008 // terms of the Library GNU General Public License as published by the
00009 // Free Software Foundation; either version 2, or (at your option)
00010 // any later version.
00011 
00012 // This library is distributed in the hope that it will be useful,
00013 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00014 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00015 // Library GNU General Public License for more details.
00016 
00017 // As a special exception, you may use this file as part of a free software
00018 // library without restriction.  Specifically, if other files instantiate
00019 // templates or use macros or inline functions from this file, or you compile
00020 // this file and link it with other files to produce an executable, this
00021 // file does not by itself cause the resulting executable to be covered by
00022 // the Library GNU General Public License.  This exception does not however
00023 // invalidate any other reasons why the executable file might be covered by
00024 // the Library GNU General Public License.
00025 
00028 #ifndef _PROINTERVAL_H_
00029 #define _PROINTERVAL_H_
00030 
00031 #include <iostream>
00032 #include <counted_ptr.h>
00033 
00034 namespace coco {
00035 
00036 class proj_rational
00037 {
00038 private:
00039   int nom;
00040   unsigned int denom;
00041 
00042 private:
00043   unsigned int abs(int n) const { return n < 0 ? -n : n; }
00044   unsigned int scm(unsigned int& n, unsigned int& m) const;
00045   unsigned int gcd(int n, unsigned int m) const;
00046   void cancel();
00047 
00048 public:
00049   proj_rational(int n=0, unsigned int d=1) : nom(n), denom(d) { if(d!=1) cancel(); }
00050 
00051   proj_rational(const proj_rational& R) : nom(R.nom), denom(R.denom) {}
00052 
00053   ~proj_rational() {}
00054 
00055   proj_rational& operator=(const proj_rational& R)
00056         { nom=R.nom; denom=R.denom; return *this; }
00057 
00058   proj_rational& operator=(int n)
00059         { nom=n; denom=1; return *this; }
00060 
00061   proj_rational operator-() const
00062         { proj_rational r(-nom,denom); return r; }
00063 
00064   proj_rational operator+(const proj_rational& R) const;
00065   proj_rational operator-(const proj_rational& R) const;
00066   proj_rational operator*(const proj_rational& R) const;
00067   proj_rational operator/(const proj_rational& R) const;
00068 
00069   proj_rational& operator+=(const proj_rational& R);
00070   proj_rational& operator-=(const proj_rational& R);
00071   proj_rational& operator*=(const proj_rational& R);
00072   proj_rational& operator/=(const proj_rational& R);
00073 
00074   interval ival() const { interval n((double)nom); return n/(double)denom; }
00075 
00076   bool operator==(const proj_rational& R) const
00077         { return nom == R.nom && denom == R.denom; }
00078 
00079   bool operator!=(const proj_rational& R) const
00080         { return nom != R.nom || denom != R.denom; }
00081 
00082   bool operator<(const proj_rational& R) const
00083         { return nom*(int)R.denom < (int)denom*R.nom; }
00084 
00085   bool operator>(const proj_rational& R) const
00086         { return nom*(int)R.denom > (int)denom*R.nom; }
00087 
00088   bool operator<=(const proj_rational& R) const
00089         { return nom*(int)R.denom <= (int)denom*R.nom; }
00090 
00091   bool operator>=(const proj_rational& R) const
00092         { return nom*(int)R.denom >= (int)denom*R.nom; }
00093 
00094   bool operator==(int n) const { return nom == n && denom == 1; }
00095 
00096   bool operator!=(int n) const { return nom != n || denom != 1; }
00097 
00098   bool operator<(int n) const { return nom < n*(int)denom; }
00099 
00100   bool operator>(int n) const { return nom > n*(int)denom; }
00101 
00102   bool operator<=(int n) const { return nom <= n*(int)denom; }
00103 
00104   bool operator>=(int n) const { return nom >= n*(int)denom; }
00105 
00106   int nominator() const { return nom; }
00107   int denominator() const { return denom; }
00108   proj_rational& set(int n, unsigned int d)
00109         { nom=n; denom=d; cancel(); return *this; }
00110 
00111   friend proj_rational abs(const proj_rational& r);
00112 };
00113 
00114 bool operator<(int n, const proj_rational& R) { return R>n; }
00115 bool operator<=(int n, const proj_rational& R) { return R>=n; }
00116 bool operator>(int n, const proj_rational& R) { return R<n; }
00117 bool operator>=(int n, const proj_rational& R) { return R<=n; }
00118 bool operator==(int n, const proj_rational& R) { return R==n; }
00119 bool operator!=(int n, const proj_rational& R) { return R!=n; }
00120 proj_rational abs(const proj_rational& r)
00121                                 { return proj_rational(r.abs(r.nom),r.denom); }
00122 
00123 template <class I>
00124 class projective_interval
00125 {
00126 private:
00127   typedef       I                      interval;
00128   typedef       projective_interval<I> proj_interval;
00129 
00130 private:
00131   interval                   i;
00132   counted_ptr<interval>      t;
00133   proj_rational              r;
00134 
00135 private:
00136   void shift(const proj_rational& s);
00137   interval shiftr(const proj_rational& s) const;
00138   interval shiftr(const proj_interval& x, const proj_rational& s) const;
00139   interval shiftintv(const proj_rational& r, const proj_rational& s) const;
00140 
00141 public:
00143   projective_interval() : i(), t(NULL), r() {}
00145   projective_interval(const proj_interval& _p) : i(_p.i), t(_p.t), r(_p.r) {}
00146 
00148   projective_interval(const interval &_t)
00149                         : i(), t(new interval(_t)), r() {}
00150   projective_interval(const interval &_i, const interval &_t)
00151                         : i(_i), t(new interval(_t)), r() {}
00152   projective_interval(const interval &_i, int _m, const interval &_t)
00153                         : i(_i), t(new interval(_t)), r(_m) {}
00154   projective_interval(const interval &_i, proj_rational _r, const interval &_t)
00155                         : i(_i), t(new interval(_t)), r(_r) {}
00156 
00158   projective_interval(const interval &_i, int _m, const proj_interval& _p)
00159                                 : i(_i), t(_p.t), r(_m) {}
00160   projective_interval(const interval &_i, proj_rational _r, const proj_interval& _p)
00161                                 : i(_i), t(_p.t), r(_r) {}
00162   projective_interval(double lo, double up, const proj_interval& _p)
00163                                 : i(lo, up), t(_p.t), r() {}
00164   projective_interval(double d, const proj_interval& _p)
00165                                 : i(d), t(_p.t), r() {}
00166   projective_interval(int d, const proj_interval& _p)
00167                                 : i((double) d), t(_p.t), r() {}
00168   projective_interval(unsigned d, const proj_interval& _p)
00169                                 : i((double) d, t(_p.t)), r() {}
00170   projective_interval(long d, const proj_interval& _p)
00171                                 : i((double) d, t(_p.t)), r() {}
00172   projective_interval(unsigned long d, const proj_interval& _p)
00173                                 : i((double) d, t(_p.t)), r() {}
00174   projective_interval(const interval& x, const proj_interval& _p)
00175                                 : i(x, t(_p.t)), r() {}
00176 
00178   virtual ~projective_interval() {}
00179 
00180 public:
00181   
00182   void set(double lo, double up, const proj_interval& _p);
00183   void set_i(double lo, double up);
00184   void set_i(const interval& j);
00185   void set_t(double lo, double up);
00186   void set_t(const interval& u);
00187   void set_r(int n, unsigned int d);
00188   void set_r(const proj_rational& s);
00189 
00190   const interval& pi() const {return i;}
00191   const interval& pt() const {return *t;}
00192   const proj_rational& pr() const {return r;}
00193   
00194   double i_inf() const {return i.inf();}
00195   double i_sup() const {return i.sup();}
00196   double t_inf() const {return t->inf();}
00197   double t_sup() const {return t->sup();}
00198   double t_exp() const {return ((double)nom)/denom;}
00199   double t_exp_n() const {return nom;}
00200   double t_exp_d() const {return denom;}
00201   
00202   static proj_interval EMPTY(const proj_interval& p)
00203         {return proj_interval(interval::empty(), interval(0,1), p);}
00204   
00205   bool empty() const {return i.empty() || t->empty();}
00206   bool is_empty() const {return i.empty() || t->empty();}
00207   bool is_i_thin() const {return i.is_thin();}
00208   bool is_t_thin() const {return t->is_thin();}
00209   bool is_unbounded_below() const
00210         {return i.inf() == -COCO_INF || (i.inf() < 0 && t->inf() == 0 && nom > 0);}
00211   bool is_unbounded_above() const
00212         {return i.sup() == -COCO_INF || (i.sup() > 0 && t->inf() == 0 && nom > 0);}
00213   bool is_entire() const {return is_unbounded_below() && is_unbounded_above();}
00214   bool is_bounded() const
00215         {return !is_unbounded_below() && !is_unbounded_above();}
00216   bool i_subset (const proj_interval& x) const {return i.subset(x.i);}
00217   bool t_subset (const proj_interval& x) const {return t->subset(*x.t);}
00218   bool certain_subset (const proj_interval& x) const;
00219   bool possible_subset (const proj_interval& x) const;   // possible subset
00220   bool i_superset (const proj_interval& x) const 
00221         {return i.superset(x.i);}
00222   bool t_superset (const proj_interval& x) const 
00223         {return t->superset(*x.t);}
00224   bool certain_superset (const proj_interval& x) const    // certain superset
00225         {return x.certain_subset(*this);}
00226   bool possible_superset (const proj_interval& x) const   // possible superset
00227         {return x.possible_subset(*this);}
00228   bool i_proper_subset (const interval& x) const 
00229         {return i.proper_subset(x.i);}
00230   bool t_proper_subset (const interval& x) const 
00231         {return t->proper_subset(*x.t);}
00232   bool certain_proper_subset (const proj_interval& x) const; //certain proper subset
00233   bool possible_proper_subset (const proj_interval& x) const;//possible proper subset
00234   bool i_proper_superset (const proj_interval& x) const 
00235         {return i.proper_superset(x.i);}
00236   bool t_proper_superset (const proj_interval& x) const 
00237         {return t->proper_superset(*x.t);}
00238   bool certain_proper_superset (const proj_interval& x) const // certain proper superset
00239         {return x.certain_proper_subset(*this);}
00240   bool possible_proper_superset (const proj_interval& x) const // possible proper superset
00241         {return x.possible_proper_subset(*this);}
00242   bool disjoint (const proj_interval& x) const;
00243 
00244   double rel_width() const { return i.rel_width(); }
00245 
00246   // y.intersect(x) intersects y with x
00247   proj_interval& intersectwith(const proj_interval& x);
00248   
00249   // y.hull(x) builds the hull of y with x
00250   proj_interval& hull(const proj_interval& x);
00251 
00252   proj_interval& hulldiff(const proj_interval& x);
00253 
00254   bool certainly_contains(const interval& x) const
00255   {
00256     return certainly_superset(proj_interval(x));
00257   }
00258 
00259   bool probably_contains(const interval& x) const
00260   {
00261     return probably_superset(proj_interval(x));
00262   }
00263 
00264   void set_bounds(double a, double b) {i.set_bounds(a, b); r=1;}
00265   void set_lb(double d) {i.set_bounds(d, i.sup()); r=1;}
00266   void set_ub(double d) {i.set_bounds(i.inf(), d); r=1;}
00267   
00268   proj_interval& operator=(const proj_interval& x){i=x.i; r=x.r; t=x.t; return *this;}
00269   proj_interval& operator=(double d) {i=d; r=1; return *this;}
00270   proj_interval& operator=(int d) {i=d; r=1; return *this;}
00271   proj_interval& operator=(unsigned d) {i=d; r=1; return *this;}
00272   proj_interval& operator=(long d) {i=d; r=1; return *this;}
00273   proj_interval& operator=(unsigned long d) {i=d; r=1; return *this;}
00274   //interval& operator=(const __Base& _x)
00275   //        {I.assign(_x.lower(),_x.upper()); return *this;}
00276 
00277   proj_interval operator-() const {return i=-i;}
00278 
00280   proj_interval& ipow(int n);
00282   proj_interval& imin(const proj_interval& x);
00284   proj_interval& imax(const proj_interval& x);
00285 
00287   proj_interval& round_to_integer();
00288   
00289   proj_interval& intersect_div(const proj_interval& __i, const proj_interval& __j);
00290   proj_interval& intersect_power(const proj_interval& __i, int __n);
00291   proj_interval& intersect_invsqr_wc(const proj_interval& __i, double __l, double __d);
00292   proj_interval& intersect_invpower_wc(const proj_interval& __i, double __l, int __n);
00293   proj_interval& intersect_invsin_wc(const proj_interval& __i, double __l, double __d);
00294   proj_interval& intersect_invcos_wc(const proj_interval& __i, double __l, double __d);
00295   proj_interval& intersect_invgauss_wc(const proj_interval& __i, double __l, double __m, double __s);
00296   
00297   double project(double __d) const;
00298 
00299   //unary +=
00300   proj_interval& operator+=(const proj_interval& a);
00301 
00302   //unary += for number
00303   proj_interval& operator+=(double a);
00304 
00305   //unary -=
00306   proj_interval& operator-=(const proj_interval& a);
00307 
00308   //unary -= for number
00309   proj_interval& operator-=(double a);
00310 
00311   //unary operator*=
00312   proj_interval& operator*=(const proj_interval& a);
00313 
00314   //unary operator*= for point
00315   proj_interval& operator*=(double a);
00316 
00330   double mid() const;
00331 
00332   /*Returns an upper bound for the diameter (width) of the interval of proj_interval, i.e.
00333 
00334   a.diam() == a.i.sup()-a.i.inf()
00335 
00336   Special cases in the extended system:
00337   <UL>
00338   <LI> a.diam() == NaN for a == [ EMPTY ]
00339   <LI> a.diam() == +INF for any infinite proj_interval
00340   </UL>
00341   */
00342   
00343   double diam() const;
00344   double width() const;
00345 
00361   double relDiam() const;
00362 
00374   double rad() const;
00375 
00386   double mig() const;
00387 
00399   double mag() const;
00400 
00412   double dist(const proj_interval& x) const;
00413 
00417   double safeguarded_mid() const;
00418 
00419   template <class J>
00420   friend projective_interval<J> operator+(const projective_interval<J>& a, const projective_interval<J>& b);
00421   template <class J>
00422   friend projective_interval<J> operator+(const projective_interval<J>& a, double b);
00423   template <class J>
00424   friend projective_interval<J> operator+(double b, const projective_interval<J>& a);
00425   template <class J>
00426   friend projective_interval<J> operator-(const projective_interval<J>& a, const projective_interval<J>& b);
00427   template <class J>
00428   friend projective_interval<J> operator-(const projective_interval<J>& a, double b);
00429   template <class J>
00430   friend projective_interval<J> operator-(double b, const projective_interval<J>& a);
00431   template <class J>
00432   friend projective_interval<J> cancel(const projective_interval<J>& a, const projective_interval<J>& b);
00433   template <class J>
00434   friend projective_interval<J> operator*(const projective_interval<J>& a, const projective_interval<J>& b);
00435   template <class J>
00436   friend projective_interval<J> operator*(const projective_interval<J>& a, double b);
00437   template <class J>
00438   friend projective_interval<J> operator*(double b, const projective_interval<J>& a);
00439   template <class J>
00440   friend projective_interval<J> operator/(const projective_interval<J>& a, const projective_interval<J>& b);
00441   template <class J>
00442   friend projective_interval<J> operator/(const projective_interval<J>& a, double b);
00443   template <class J>
00444   friend projective_interval<J> operator/(double b, const projective_interval<J>& a);
00445   template <class J>
00446   friend projective_interval<J> acos(const projective_interval<J>& x);
00447   template <class J>
00448   friend projective_interval<J> abs(const projective_interval<J>& x);
00449   template <class J>
00450   friend projective_interval<J> acosh(const projective_interval<J>& x);
00451   template <class J>
00452   friend projective_interval<J> acot(const projective_interval<J>& x);
00453   template <class J>
00454   friend projective_interval<J> acoth(const projective_interval<J>& x);
00455   template <class J>
00456   friend projective_interval<J> asin(const projective_interval<J>& x);
00457   template <class J>
00458   friend projective_interval<J> asinh(const projective_interval<J>& x);
00459   template <class J>
00460   friend projective_interval<J> atan(const projective_interval<J>& x);
00461   template <class J>
00462   friend projective_interval<J> atanh(const projective_interval<J>& x);
00463   template <class J>
00464   friend projective_interval<J> cos(const projective_interval<J>& x);
00465   template <class J>
00466   friend projective_interval<J> cosh(const projective_interval<J>& x);
00467   template <class J>
00468   friend projective_interval<J> cot(const projective_interval<J>& x);
00469   template <class J>
00470   friend projective_interval<J> coth(const projective_interval<J>& x);
00471   template <class J>
00472   friend projective_interval<J> exp(const projective_interval<J>& x);
00473   template <class J>
00474   friend projective_interval<J> exp10(const projective_interval<J>& x);
00475   template <class J>
00476   friend projective_interval<J> exp2(const projective_interval<J>& x);
00477   template <class J>
00478   friend projective_interval<J> expm1(const projective_interval<J>& x);
00479   template <class J>
00480   friend projective_interval<J> log(const projective_interval<J>& x);
00481   template <class J>
00482   friend projective_interval<J> log10(const projective_interval<J>& x);
00483   template <class J>
00484   friend projective_interval<J> log1p(const projective_interval<J>& x);
00485   template <class J>
00486   friend projective_interval<J> log2(const projective_interval<J>& x);
00487   template <class J>
00488   friend projective_interval<J> power(const projective_interval<J>& x, int n);
00489   template <class J>
00490   friend projective_interval<J> pow(const projective_interval<J>& x, const projective_interval<J>& y);
00491   template <class J>
00492   friend projective_interval<J> sin(const projective_interval<J>& x);
00493   template <class J>
00494   friend projective_interval<J> sinh(const projective_interval<J>& x);
00495   template <class J>
00496   friend projective_interval<J> sqr(const projective_interval<J>& x);
00497   template <class J>
00498   friend projective_interval<J> sqrt(const projective_interval<J>& x);
00499   template <class J>
00500   friend projective_interval<J> tan(const projective_interval<J>& x);
00501   template <class J>
00502   friend projective_interval<J> tanh(const projective_interval<J>& x);
00503   template <class J>
00504   friend projective_interval<J> imax(const projective_interval<J>& x, const projective_interval<J>& y);
00505   template <class J>
00506   friend projective_interval<J> imin(const projective_interval<J>& x, const projective_interval<J>& y);
00507   template <class J>
00508   friend projective_interval<J> division_part1(const projective_interval<J>& x,
00509         const projective_interval<J>& y, bool& b);
00510   template <class J>
00511   friend projective_interval<J> division_part2(const projective_interval<J>& x,
00512         const projective_interval<J>& y, bool b);
00513 };
00514 
00515 } // namespace coco
00516 
00517 #include <prointerval.hpp>
00518 
00519 #endif // _PROINTERVAL_H_