lemon/dim2.h
author Peter Kovacs <kpeter@inf.elte.hu>
Sat, 20 Feb 2010 18:39:03 +0100
changeset 839 f3bc4e9b5f3a
parent 440 88ed40ad0d4f
permissions -rw-r--r--
New heuristics for MCF algorithms (#340)
and some implementation improvements.

- A useful heuristic is added to NetworkSimplex to make the
initial pivots faster.
- A powerful global update heuristic is added to CostScaling
and the implementation is reworked with various improvements.
- Better relabeling in CostScaling to improve numerical stability
and make the code faster.
- A small improvement is made in CapacityScaling for better
delta computation.
- Add notes to the classes about the usage of vector<char> instead
of vector<bool> for efficiency reasons.
     1 /* -*- mode: C++; indent-tabs-mode: nil; -*-
     2  *
     3  * This file is a part of LEMON, a generic C++ optimization library.
     4  *
     5  * Copyright (C) 2003-2009
     6  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     7  * (Egervary Research Group on Combinatorial Optimization, EGRES).
     8  *
     9  * Permission to use, modify and distribute this software is granted
    10  * provided that this copyright notice appears in all copies. For
    11  * precise terms see the accompanying LICENSE file.
    12  *
    13  * This software is provided "AS IS" with no warranty of any kind,
    14  * express or implied, and with no claim as to its suitability for any
    15  * purpose.
    16  *
    17  */
    18 
    19 #ifndef LEMON_DIM2_H
    20 #define LEMON_DIM2_H
    21 
    22 #include <iostream>
    23 
    24 ///\ingroup geomdat
    25 ///\file
    26 ///\brief A simple two dimensional vector and a bounding box implementation
    27 
    28 namespace lemon {
    29 
    30   ///Tools for handling two dimensional coordinates
    31 
    32   ///This namespace is a storage of several
    33   ///tools for handling two dimensional coordinates
    34   namespace dim2 {
    35 
    36   /// \addtogroup geomdat
    37   /// @{
    38 
    39   /// Two dimensional vector (plain vector)
    40 
    41   /// A simple two dimensional vector (plain vector) implementation
    42   /// with the usual vector operations.
    43   template<typename T>
    44     class Point {
    45 
    46     public:
    47 
    48       typedef T Value;
    49 
    50       ///First coordinate
    51       T x;
    52       ///Second coordinate
    53       T y;
    54 
    55       ///Default constructor
    56       Point() {}
    57 
    58       ///Construct an instance from coordinates
    59       Point(T a, T b) : x(a), y(b) { }
    60 
    61       ///Returns the dimension of the vector (i.e. returns 2).
    62 
    63       ///The dimension of the vector.
    64       ///This function always returns 2.
    65       int size() const { return 2; }
    66 
    67       ///Subscripting operator
    68 
    69       ///\c p[0] is \c p.x and \c p[1] is \c p.y
    70       ///
    71       T& operator[](int idx) { return idx == 0 ? x : y; }
    72 
    73       ///Const subscripting operator
    74 
    75       ///\c p[0] is \c p.x and \c p[1] is \c p.y
    76       ///
    77       const T& operator[](int idx) const { return idx == 0 ? x : y; }
    78 
    79       ///Conversion constructor
    80       template<class TT> Point(const Point<TT> &p) : x(p.x), y(p.y) {}
    81 
    82       ///Give back the square of the norm of the vector
    83       T normSquare() const {
    84         return x*x+y*y;
    85       }
    86 
    87       ///Increment the left hand side by \c u
    88       Point<T>& operator +=(const Point<T>& u) {
    89         x += u.x;
    90         y += u.y;
    91         return *this;
    92       }
    93 
    94       ///Decrement the left hand side by \c u
    95       Point<T>& operator -=(const Point<T>& u) {
    96         x -= u.x;
    97         y -= u.y;
    98         return *this;
    99       }
   100 
   101       ///Multiply the left hand side with a scalar
   102       Point<T>& operator *=(const T &u) {
   103         x *= u;
   104         y *= u;
   105         return *this;
   106       }
   107 
   108       ///Divide the left hand side by a scalar
   109       Point<T>& operator /=(const T &u) {
   110         x /= u;
   111         y /= u;
   112         return *this;
   113       }
   114 
   115       ///Return the scalar product of two vectors
   116       T operator *(const Point<T>& u) const {
   117         return x*u.x+y*u.y;
   118       }
   119 
   120       ///Return the sum of two vectors
   121       Point<T> operator+(const Point<T> &u) const {
   122         Point<T> b=*this;
   123         return b+=u;
   124       }
   125 
   126       ///Return the negative of the vector
   127       Point<T> operator-() const {
   128         Point<T> b=*this;
   129         b.x=-b.x; b.y=-b.y;
   130         return b;
   131       }
   132 
   133       ///Return the difference of two vectors
   134       Point<T> operator-(const Point<T> &u) const {
   135         Point<T> b=*this;
   136         return b-=u;
   137       }
   138 
   139       ///Return a vector multiplied by a scalar
   140       Point<T> operator*(const T &u) const {
   141         Point<T> b=*this;
   142         return b*=u;
   143       }
   144 
   145       ///Return a vector divided by a scalar
   146       Point<T> operator/(const T &u) const {
   147         Point<T> b=*this;
   148         return b/=u;
   149       }
   150 
   151       ///Test equality
   152       bool operator==(const Point<T> &u) const {
   153         return (x==u.x) && (y==u.y);
   154       }
   155 
   156       ///Test inequality
   157       bool operator!=(Point u) const {
   158         return  (x!=u.x) || (y!=u.y);
   159       }
   160 
   161     };
   162 
   163   ///Return a Point
   164 
   165   ///Return a Point.
   166   ///\relates Point
   167   template <typename T>
   168   inline Point<T> makePoint(const T& x, const T& y) {
   169     return Point<T>(x, y);
   170   }
   171 
   172   ///Return a vector multiplied by a scalar
   173 
   174   ///Return a vector multiplied by a scalar.
   175   ///\relates Point
   176   template<typename T> Point<T> operator*(const T &u,const Point<T> &x) {
   177     return x*u;
   178   }
   179 
   180   ///Read a plain vector from a stream
   181 
   182   ///Read a plain vector from a stream.
   183   ///\relates Point
   184   ///
   185   template<typename T>
   186   inline std::istream& operator>>(std::istream &is, Point<T> &z) {
   187     char c;
   188     if (is >> c) {
   189       if (c != '(') is.putback(c);
   190     } else {
   191       is.clear();
   192     }
   193     if (!(is >> z.x)) return is;
   194     if (is >> c) {
   195       if (c != ',') is.putback(c);
   196     } else {
   197       is.clear();
   198     }
   199     if (!(is >> z.y)) return is;
   200     if (is >> c) {
   201       if (c != ')') is.putback(c);
   202     } else {
   203       is.clear();
   204     }
   205     return is;
   206   }
   207 
   208   ///Write a plain vector to a stream
   209 
   210   ///Write a plain vector to a stream.
   211   ///\relates Point
   212   ///
   213   template<typename T>
   214   inline std::ostream& operator<<(std::ostream &os, const Point<T>& z)
   215   {
   216     os << "(" << z.x << "," << z.y << ")";
   217     return os;
   218   }
   219 
   220   ///Rotate by 90 degrees
   221 
   222   ///Returns the parameter rotated by 90 degrees in positive direction.
   223   ///\relates Point
   224   ///
   225   template<typename T>
   226   inline Point<T> rot90(const Point<T> &z)
   227   {
   228     return Point<T>(-z.y,z.x);
   229   }
   230 
   231   ///Rotate by 180 degrees
   232 
   233   ///Returns the parameter rotated by 180 degrees.
   234   ///\relates Point
   235   ///
   236   template<typename T>
   237   inline Point<T> rot180(const Point<T> &z)
   238   {
   239     return Point<T>(-z.x,-z.y);
   240   }
   241 
   242   ///Rotate by 270 degrees
   243 
   244   ///Returns the parameter rotated by 90 degrees in negative direction.
   245   ///\relates Point
   246   ///
   247   template<typename T>
   248   inline Point<T> rot270(const Point<T> &z)
   249   {
   250     return Point<T>(z.y,-z.x);
   251   }
   252 
   253 
   254 
   255   /// Bounding box of plain vectors (points).
   256 
   257   /// A class to calculate or store the bounding box of plain vectors
   258   /// (\ref Point "points").
   259   template<typename T>
   260   class Box {
   261       Point<T> _bottom_left, _top_right;
   262       bool _empty;
   263     public:
   264 
   265       ///Default constructor: creates an empty box
   266       Box() { _empty = true; }
   267 
   268       ///Construct a box from one point
   269       Box(Point<T> a) {
   270         _bottom_left = _top_right = a;
   271         _empty = false;
   272       }
   273 
   274       ///Construct a box from two points
   275 
   276       ///Construct a box from two points.
   277       ///\param a The bottom left corner.
   278       ///\param b The top right corner.
   279       ///\warning The coordinates of the bottom left corner must be no more
   280       ///than those of the top right one.
   281       Box(Point<T> a,Point<T> b)
   282       {
   283         _bottom_left = a;
   284         _top_right = b;
   285         _empty = false;
   286       }
   287 
   288       ///Construct a box from four numbers
   289 
   290       ///Construct a box from four numbers.
   291       ///\param l The left side of the box.
   292       ///\param b The bottom of the box.
   293       ///\param r The right side of the box.
   294       ///\param t The top of the box.
   295       ///\warning The left side must be no more than the right side and
   296       ///bottom must be no more than the top.
   297       Box(T l,T b,T r,T t)
   298       {
   299         _bottom_left=Point<T>(l,b);
   300         _top_right=Point<T>(r,t);
   301         _empty = false;
   302       }
   303 
   304       ///Return \c true if the box is empty.
   305 
   306       ///Return \c true if the box is empty (i.e. return \c false
   307       ///if at least one point was added to the box or the coordinates of
   308       ///the box were set).
   309       ///
   310       ///The coordinates of an empty box are not defined.
   311       bool empty() const {
   312         return _empty;
   313       }
   314 
   315       ///Make the box empty
   316       void clear() {
   317         _empty = true;
   318       }
   319 
   320       ///Give back the bottom left corner of the box
   321 
   322       ///Give back the bottom left corner of the box.
   323       ///If the box is empty, then the return value is not defined.
   324       Point<T> bottomLeft() const {
   325         return _bottom_left;
   326       }
   327 
   328       ///Set the bottom left corner of the box
   329 
   330       ///Set the bottom left corner of the box.
   331       ///\pre The box must not be empty.
   332       void bottomLeft(Point<T> p) {
   333         _bottom_left = p;
   334       }
   335 
   336       ///Give back the top right corner of the box
   337 
   338       ///Give back the top right corner of the box.
   339       ///If the box is empty, then the return value is not defined.
   340       Point<T> topRight() const {
   341         return _top_right;
   342       }
   343 
   344       ///Set the top right corner of the box
   345 
   346       ///Set the top right corner of the box.
   347       ///\pre The box must not be empty.
   348       void topRight(Point<T> p) {
   349         _top_right = p;
   350       }
   351 
   352       ///Give back the bottom right corner of the box
   353 
   354       ///Give back the bottom right corner of the box.
   355       ///If the box is empty, then the return value is not defined.
   356       Point<T> bottomRight() const {
   357         return Point<T>(_top_right.x,_bottom_left.y);
   358       }
   359 
   360       ///Set the bottom right corner of the box
   361 
   362       ///Set the bottom right corner of the box.
   363       ///\pre The box must not be empty.
   364       void bottomRight(Point<T> p) {
   365         _top_right.x = p.x;
   366         _bottom_left.y = p.y;
   367       }
   368 
   369       ///Give back the top left corner of the box
   370 
   371       ///Give back the top left corner of the box.
   372       ///If the box is empty, then the return value is not defined.
   373       Point<T> topLeft() const {
   374         return Point<T>(_bottom_left.x,_top_right.y);
   375       }
   376 
   377       ///Set the top left corner of the box
   378 
   379       ///Set the top left corner of the box.
   380       ///\pre The box must not be empty.
   381       void topLeft(Point<T> p) {
   382         _top_right.y = p.y;
   383         _bottom_left.x = p.x;
   384       }
   385 
   386       ///Give back the bottom of the box
   387 
   388       ///Give back the bottom of the box.
   389       ///If the box is empty, then the return value is not defined.
   390       T bottom() const {
   391         return _bottom_left.y;
   392       }
   393 
   394       ///Set the bottom of the box
   395 
   396       ///Set the bottom of the box.
   397       ///\pre The box must not be empty.
   398       void bottom(T t) {
   399         _bottom_left.y = t;
   400       }
   401 
   402       ///Give back the top of the box
   403 
   404       ///Give back the top of the box.
   405       ///If the box is empty, then the return value is not defined.
   406       T top() const {
   407         return _top_right.y;
   408       }
   409 
   410       ///Set the top of the box
   411 
   412       ///Set the top of the box.
   413       ///\pre The box must not be empty.
   414       void top(T t) {
   415         _top_right.y = t;
   416       }
   417 
   418       ///Give back the left side of the box
   419 
   420       ///Give back the left side of the box.
   421       ///If the box is empty, then the return value is not defined.
   422       T left() const {
   423         return _bottom_left.x;
   424       }
   425 
   426       ///Set the left side of the box
   427 
   428       ///Set the left side of the box.
   429       ///\pre The box must not be empty.
   430       void left(T t) {
   431         _bottom_left.x = t;
   432       }
   433 
   434       /// Give back the right side of the box
   435 
   436       /// Give back the right side of the box.
   437       ///If the box is empty, then the return value is not defined.
   438       T right() const {
   439         return _top_right.x;
   440       }
   441 
   442       ///Set the right side of the box
   443 
   444       ///Set the right side of the box.
   445       ///\pre The box must not be empty.
   446       void right(T t) {
   447         _top_right.x = t;
   448       }
   449 
   450       ///Give back the height of the box
   451 
   452       ///Give back the height of the box.
   453       ///If the box is empty, then the return value is not defined.
   454       T height() const {
   455         return _top_right.y-_bottom_left.y;
   456       }
   457 
   458       ///Give back the width of the box
   459 
   460       ///Give back the width of the box.
   461       ///If the box is empty, then the return value is not defined.
   462       T width() const {
   463         return _top_right.x-_bottom_left.x;
   464       }
   465 
   466       ///Checks whether a point is inside the box
   467       bool inside(const Point<T>& u) const {
   468         if (_empty)
   469           return false;
   470         else {
   471           return ( (u.x-_bottom_left.x)*(_top_right.x-u.x) >= 0 &&
   472                    (u.y-_bottom_left.y)*(_top_right.y-u.y) >= 0 );
   473         }
   474       }
   475 
   476       ///Increments the box with a point
   477 
   478       ///Increments the box with a point.
   479       ///
   480       Box& add(const Point<T>& u){
   481         if (_empty) {
   482           _bottom_left = _top_right = u;
   483           _empty = false;
   484         }
   485         else {
   486           if (_bottom_left.x > u.x) _bottom_left.x = u.x;
   487           if (_bottom_left.y > u.y) _bottom_left.y = u.y;
   488           if (_top_right.x < u.x) _top_right.x = u.x;
   489           if (_top_right.y < u.y) _top_right.y = u.y;
   490         }
   491         return *this;
   492       }
   493 
   494       ///Increments the box to contain another box
   495 
   496       ///Increments the box to contain another box.
   497       ///
   498       Box& add(const Box &u){
   499         if ( !u.empty() ){
   500           add(u._bottom_left);
   501           add(u._top_right);
   502         }
   503         return *this;
   504       }
   505 
   506       ///Intersection of two boxes
   507 
   508       ///Intersection of two boxes.
   509       ///
   510       Box operator&(const Box& u) const {
   511         Box b;
   512         if (_empty || u._empty) {
   513           b._empty = true;
   514         } else {
   515           b._bottom_left.x = std::max(_bottom_left.x, u._bottom_left.x);
   516           b._bottom_left.y = std::max(_bottom_left.y, u._bottom_left.y);
   517           b._top_right.x = std::min(_top_right.x, u._top_right.x);
   518           b._top_right.y = std::min(_top_right.y, u._top_right.y);
   519           b._empty = b._bottom_left.x > b._top_right.x ||
   520                      b._bottom_left.y > b._top_right.y;
   521         }
   522         return b;
   523       }
   524 
   525   };//class Box
   526 
   527 
   528   ///Read a box from a stream
   529 
   530   ///Read a box from a stream.
   531   ///\relates Box
   532   template<typename T>
   533   inline std::istream& operator>>(std::istream &is, Box<T>& b) {
   534     char c;
   535     Point<T> p;
   536     if (is >> c) {
   537       if (c != '(') is.putback(c);
   538     } else {
   539       is.clear();
   540     }
   541     if (!(is >> p)) return is;
   542     b.bottomLeft(p);
   543     if (is >> c) {
   544       if (c != ',') is.putback(c);
   545     } else {
   546       is.clear();
   547     }
   548     if (!(is >> p)) return is;
   549     b.topRight(p);
   550     if (is >> c) {
   551       if (c != ')') is.putback(c);
   552     } else {
   553       is.clear();
   554     }
   555     return is;
   556   }
   557 
   558   ///Write a box to a stream
   559 
   560   ///Write a box to a stream.
   561   ///\relates Box
   562   template<typename T>
   563   inline std::ostream& operator<<(std::ostream &os, const Box<T>& b)
   564   {
   565     os << "(" << b.bottomLeft() << "," << b.topRight() << ")";
   566     return os;
   567   }
   568 
   569   ///Map of x-coordinates of a <tt>Point</tt>-map
   570 
   571   ///Map of x-coordinates of a \ref Point "Point"-map.
   572   ///
   573   template<class M>
   574   class XMap
   575   {
   576     M& _map;
   577   public:
   578 
   579     typedef typename M::Value::Value Value;
   580     typedef typename M::Key Key;
   581     ///\e
   582     XMap(M& map) : _map(map) {}
   583     Value operator[](Key k) const {return _map[k].x;}
   584     void set(Key k,Value v) {_map.set(k,typename M::Value(v,_map[k].y));}
   585   };
   586 
   587   ///Returns an XMap class
   588 
   589   ///This function just returns an XMap class.
   590   ///\relates XMap
   591   template<class M>
   592   inline XMap<M> xMap(M &m)
   593   {
   594     return XMap<M>(m);
   595   }
   596 
   597   template<class M>
   598   inline XMap<M> xMap(const M &m)
   599   {
   600     return XMap<M>(m);
   601   }
   602 
   603   ///Constant (read only) version of XMap
   604 
   605   ///Constant (read only) version of XMap.
   606   ///
   607   template<class M>
   608   class ConstXMap
   609   {
   610     const M& _map;
   611   public:
   612 
   613     typedef typename M::Value::Value Value;
   614     typedef typename M::Key Key;
   615     ///\e
   616     ConstXMap(const M &map) : _map(map) {}
   617     Value operator[](Key k) const {return _map[k].x;}
   618   };
   619 
   620   ///Returns a ConstXMap class
   621 
   622   ///This function just returns a ConstXMap class.
   623   ///\relates ConstXMap
   624   template<class M>
   625   inline ConstXMap<M> xMap(const M &m)
   626   {
   627     return ConstXMap<M>(m);
   628   }
   629 
   630   ///Map of y-coordinates of a <tt>Point</tt>-map
   631 
   632   ///Map of y-coordinates of a \ref Point "Point"-map.
   633   ///
   634   template<class M>
   635   class YMap
   636   {
   637     M& _map;
   638   public:
   639 
   640     typedef typename M::Value::Value Value;
   641     typedef typename M::Key Key;
   642     ///\e
   643     YMap(M& map) : _map(map) {}
   644     Value operator[](Key k) const {return _map[k].y;}
   645     void set(Key k,Value v) {_map.set(k,typename M::Value(_map[k].x,v));}
   646   };
   647 
   648   ///Returns a YMap class
   649 
   650   ///This function just returns a YMap class.
   651   ///\relates YMap
   652   template<class M>
   653   inline YMap<M> yMap(M &m)
   654   {
   655     return YMap<M>(m);
   656   }
   657 
   658   template<class M>
   659   inline YMap<M> yMap(const M &m)
   660   {
   661     return YMap<M>(m);
   662   }
   663 
   664   ///Constant (read only) version of YMap
   665 
   666   ///Constant (read only) version of YMap.
   667   ///
   668   template<class M>
   669   class ConstYMap
   670   {
   671     const M& _map;
   672   public:
   673 
   674     typedef typename M::Value::Value Value;
   675     typedef typename M::Key Key;
   676     ///\e
   677     ConstYMap(const M &map) : _map(map) {}
   678     Value operator[](Key k) const {return _map[k].y;}
   679   };
   680 
   681   ///Returns a ConstYMap class
   682 
   683   ///This function just returns a ConstYMap class.
   684   ///\relates ConstYMap
   685   template<class M>
   686   inline ConstYMap<M> yMap(const M &m)
   687   {
   688     return ConstYMap<M>(m);
   689   }
   690 
   691 
   692   ///\brief Map of the normSquare() of a <tt>Point</tt>-map
   693   ///
   694   ///Map of the \ref Point::normSquare() "normSquare()"
   695   ///of a \ref Point "Point"-map.
   696   template<class M>
   697   class NormSquareMap
   698   {
   699     const M& _map;
   700   public:
   701 
   702     typedef typename M::Value::Value Value;
   703     typedef typename M::Key Key;
   704     ///\e
   705     NormSquareMap(const M &map) : _map(map) {}
   706     Value operator[](Key k) const {return _map[k].normSquare();}
   707   };
   708 
   709   ///Returns a NormSquareMap class
   710 
   711   ///This function just returns a NormSquareMap class.
   712   ///\relates NormSquareMap
   713   template<class M>
   714   inline NormSquareMap<M> normSquareMap(const M &m)
   715   {
   716     return NormSquareMap<M>(m);
   717   }
   718 
   719   /// @}
   720 
   721   } //namespce dim2
   722 
   723 } //namespace lemon
   724 
   725 #endif //LEMON_DIM2_H