Fixed nasty bug in BipartiteGraph::erase1 and erase2
[libdai.git] / include / dai / bipgraph.h
index d7decba..c9701ef 100644 (file)
 #define __defined_libdai_bipgraph_h
 
 
+#include <ostream>
 #include <vector>
+#include <cassert>
 #include <algorithm>
-#include <boost/numeric/ublas/vector_sparse.hpp>
+#include <dai/util.h>
 
 
 namespace dai {
 
 
-/// A BipartiteGraph represents a graph with two types of nodes
-template <class T1, class T2> class BipartiteGraph {
+/// A BipartiteGraph represents the neighborhood structure of nodes in a bipartite graph.
+/** A bipartite graph has two types of nodes: type 1 and type 2. Edges can occur only between
+ *  nodes of different type. Nodes are indexed by an unsigned integer, edges are indexed as
+ *  a pair of unsigned integers (where the pair (a,b) means the b'th neighbor of the a'th node).
+ *  The BipartiteGraph stores neighborhood structures as vectors of vectors of Neighbor entries.
+ */
+class BipartiteGraph {
     public:
-        typedef std::pair<size_t,size_t>    _edge_t;
-        typedef std::vector<size_t>         _nb_t;
-        typedef _nb_t::const_iterator       _nb_cit;
-        
+        /// A Neighbor describes a neighboring node of some other node.
+        /** Iterating over all neighbors of node n1 of type 1 can be done in the following way:
+         *  \code
+         *      foreach( const BipartiteGraph::Neighbor &n2, nb1(n1) ) {
+         *          size_t _n2 = n2.iter;
+         *          size_t _n1 = n2.dual;
+         *          // n2 == n2.node;
+         *          // The _n2'th neighbor of n1 is n2, and the _n1'th neighbor of n2 is n1:
+         *          // nb1(n1)[_n2] == n2, nb2(n2)[_n1] == n1 
+         *      }
+         *  \endcode
+         */
+        struct Neighbor {
+            /// Corresponds to the index of this Neighbor entry in the vector of neighbors
+            unsigned iter;
+            /// Contains the number of the neighboring node
+            unsigned node;
+            /// Contains the "dual" iter
+            unsigned dual;
+            /// Cast to unsigned returns node
+            operator unsigned () const { return node; }
+            /// Default constructor
+            Neighbor() {}
+            /// Constructor
+            Neighbor( size_t iter, size_t node, size_t dual ) : iter(iter), node(node), dual(dual) {}
+        };
+
+        /// Neighbors is a vector of Neighbor entries; each node has an associated Neighbors variable, which describes its neighbors.
+        typedef std::vector<Neighbor> Neighbors;
+
+        /// Edge is used as index of an edge: an Edge(a,b) corresponds to the edge between the a'th node and its b'th neighbor.
+        typedef std::pair<size_t,size_t> Edge;
 
     private:
-        /// Vertices of first type
-        std::vector<T1>                     _V1;
-
-        /// Vertices of second type
-        std::vector<T2>                     _V2;
-        
-        /// Edges, which are pairs (v1,v2) with v1 in _V1 and v2 in _V2
-        std::vector<_edge_t>                _E12;
-
-        /// Conversion matrix that computes which index of _E12 corresponds to a vertex index pair (v1,v2)
-        std::vector<boost::numeric::ublas::mapped_vector<size_t> > _E12ind;
-
-        /// Neighbour indices of vertices of first type
-        std::vector<_nb_t>                  _nb1;
-
-        /// Neighbour indices of vertices of second type
-        std::vector<_nb_t>                  _nb2;
+        /// _nb1 contains for each node of type 1 a vector of its neighbors
+        std::vector<Neighbors> _nb1;
+        /// _nb2 contains for each node of type 2 a vector of its neighbors
+        std::vector<Neighbors> _nb2;
 
+        /// Used internally by isTree()
+        struct levelType {
+            std::vector<size_t> ind1;       // indices of vertices of type 1
+            std::vector<size_t> ind2;       // indices of vertices of type 2
+        };
 
     public:
         /// Default constructor
-        BipartiteGraph<T1,T2> () {};
+        BipartiteGraph() : _nb1(), _nb2() {}
 
         /// Copy constructor
-        BipartiteGraph<T1,T2> ( const BipartiteGraph<T1,T2> & x ) : _V1(x._V1), _V2(x._V2), _E12(x._E12), _E12ind(x._E12ind), _nb1(x._nb1), _nb2(x._nb2) {};
+        BipartiteGraph( const BipartiteGraph & x ) : _nb1(x._nb1), _nb2(x._nb2) {}
 
         /// Assignment operator
-        BipartiteGraph<T1,T2> & operator=(const BipartiteGraph<T1,T2> & x) {
+        BipartiteGraph & operator=( const BipartiteGraph & x ) {
             if( this != &x ) {
-                _V1 =       x._V1;
-                _V2 =       x._V2;
-                _E12 =      x._E12;
-                _E12ind =   x._E12ind;
-                _nb1 =      x._nb1;
-                _nb2 =      x._nb2;
+                _nb1 = x._nb1;
+                _nb2 = x._nb2;
             }
             return *this;
         }
+
+        /// Create bipartite graph from a range of edges. 
+        /** nr1 is the number of nodes of type 1, nr2 the number of nodes of type 2. 
+         *  The value_type of an EdgeInputIterator should be Edge.
+         */
+        template<typename EdgeInputIterator>
+        void create( size_t nr1, size_t nr2, EdgeInputIterator begin, EdgeInputIterator end );
+
+        /// Construct bipartite graph from a range of edges. 
+        /** nr1 is the number of nodes of type 1, nr2 the number of nodes of type 2. 
+         *  The value_type of an EdgeInputIterator should be Edge.
+         */
+        template<typename EdgeInputIterator>
+        BipartiteGraph( size_t nr1, size_t nr2, EdgeInputIterator begin, EdgeInputIterator end ) : _nb1( nr1 ), _nb2( nr2 ) {
+            create( nr1, nr2, begin, end );
+        }
+
+        /// Returns constant reference to the _i2'th neighbor of node i1 of type 1
+        const Neighbor & nb1( size_t i1, size_t _i2 ) const { 
+#ifdef DAI_DEBUG
+            assert( i1 < _nb1.size() );
+            assert( _i2 < _nb1[i1].size() );
+#endif
+            return _nb1[i1][_i2]; 
+        }
+        /// Returns reference to the _i2'th neighbor of node i1 of type 1
+        Neighbor & nb1( size_t i1, size_t _i2 ) {
+#ifdef DAI_DEBUG
+            assert( i1 < _nb1.size() );
+            assert( _i2 < _nb1[i1].size() );
+#endif
+            return _nb1[i1][_i2]; 
+        }
+
+        /// Returns constant reference to the _i1'th neighbor of node i2 of type 2
+        const Neighbor & nb2( size_t i2, size_t _i1 ) const { 
+#ifdef DAI_DEBUG
+            assert( i2 < _nb2.size() );
+            assert( _i1 < _nb2[i2].size() );
+#endif
+            return _nb2[i2][_i1]; 
+        }
+        /// Returns reference to the _i1'th neighbor of node i2 of type 2
+        Neighbor & nb2( size_t i2, size_t _i1 ) { 
+#ifdef DAI_DEBUG
+            assert( i2 < _nb2.size() );
+            assert( _i1 < _nb2[i2].size() );
+#endif
+            return _nb2[i2][_i1]; 
+        }
+
+        /// Returns constant reference to all neighbors of node i1 of type 1
+        const Neighbors & nb1( size_t i1 ) const { 
+#ifdef DAI_DEBUG
+            assert( i1 < _nb1.size() );
+#endif
+            return _nb1[i1]; 
+        }
+        /// Returns reference to all neighbors of node of i1 type 1
+        Neighbors & nb1( size_t i1 ) { 
+#ifdef DAI_DEBUG
+            assert( i1 < _nb1.size() );
+#endif
+            return _nb1[i1]; 
+        }
+
+        /// Returns constant reference to all neighbors of node i2 of type 2
+        const Neighbors & nb2( size_t i2 ) const { 
+#ifdef DAI_DEBUG
+            assert( i2 < _nb2.size() );
+#endif
+            return _nb2[i2]; 
+        }
+        /// Returns reference to all neighbors of node i2 of type 2
+        Neighbors & nb2( size_t i2 ) { 
+#ifdef DAI_DEBUG
+            assert( i2 < _nb2.size() );
+#endif
+            return _nb2[i2]; 
+        }
+
+        /// Returns number of nodes of type 1
+        size_t nr1() const { return _nb1.size(); }
+        /// Returns number of nodes of type 2
+        size_t nr2() const { return _nb2.size(); }
         
-        /// Provides read access to node of first type
-        const T1 & V1( size_t i1 ) const { return _V1[i1]; }
-        /// Provides full access to node of first type
-        T1 & V1( size_t i1 ) { return _V1[i1]; }
-        /// Provides read access to all nodes of first type
-        const std::vector<T1> & V1s() const { return _V1; }
-        /// Provides full access to all nodes of first type
-        std::vector<T1> & V1s() { return _V1; }
-
-        /// Provides read access to node of second type
-        const T2 & V2( size_t i2 ) const { return _V2[i2]; }
-        /// Provides full access to node of second type
-        T2 & V2( size_t i2 ) { return _V2[i2]; }
-        /// Provides read access to all nodes of second type
-        const std::vector<T2> & V2s() const { return _V2; }
-        /// Provides full access to all nodes of second type
-        std::vector<T2> & V2s() { return _V2; }
-
-        /// Provides read access to edge
-        const _edge_t & edge(size_t ind) const { return _E12[ind]; }
-        /// Provides full access to edge
-        _edge_t & edge(size_t ind) { return _E12[ind]; }
-        /// Provides read access to all edges
-        const std::vector<_edge_t> & edges() const { return _E12; }
-        /// Provides full access to all edges
-        std::vector<_edge_t> & edges() { return _E12; }
-        /// Returns number of edges
-        size_t nr_edges() const { return _E12.size(); }
-
-        /// Provides read access to neighbours of node of first type
-        const _nb_t & nb1( size_t i1 ) const { return _nb1[i1]; }
-        /// Provides full access to neighbours of node of first type
-        _nb_t & nb1( size_t i1 ) { return _nb1[i1]; }
-
-        /// Provides read access to neighbours of node of second type
-        const _nb_t & nb2( size_t i2 ) const { return _nb2[i2]; }
-        /// Provides full access to neighbours of node of second type
-        _nb_t & nb2( size_t i2 ) { return _nb2[i2]; }
-
-        /// Converts the pair of indices (i1,i2) to the corresponding edge index
-        size_t VV2E( const size_t i1, const size_t i2 ) const { return _E12ind[i1](i2); }
-
-
-        /// Regenerates internal structures (_E12ind, _nb1 and _nb2) based upon _V1, _V2 and _E12
-        void Regenerate() {
-            // Calculate _nb1 and _nb2
-
-            // Start with empty vectors
-            _nb1.clear();
-            _nb1.resize(_V1.size());
-            // Start with empty vectors
-            _nb2.clear();
-            _nb2.resize(_V2.size());
-            // Each edge yields a neighbour pair
-            for( std::vector<_edge_t>::const_iterator e = _E12.begin(); e != _E12.end(); e++ ) {
-                _nb1[e->first].push_back(e->second);
-                _nb2[e->second].push_back(e->first);
-            }
-            // Remove duplicates from _nb1
-            for( size_t i1 = 0; i1 < _V1.size(); i1++ ) {
-                _nb_t::iterator new_end = unique(_nb1[i1].begin(), _nb1[i1].end());
-                _nb1[i1].erase( new_end, _nb1[i1].end() );
-            }
-            // Remove duplicates from _nb2
-            for( size_t i2 = 0; i2 < _V2.size(); i2++ ) {
-                _nb_t::iterator new_end = unique(_nb2[i2].begin(), _nb2[i2].end());
-                _nb2[i2].erase( new_end, _nb2[i2].end() );
+        /// Calculates the number of edges
+        size_t nrEdges() const {
+            size_t sum = 0;
+            for( size_t i1 = 0; i1 < nr1(); i1++ )
+                sum += nb1(i1).size();
+            return sum;
+        }
+        
+        /// Add node of type 1 without neighbors.
+        void add1() { _nb1.push_back( Neighbors() ); }
+        
+        /// Add node of type 2 without neighbors.
+        void add2() { _nb2.push_back( Neighbors() ); }
+
+        /// Add node of type 1 with neighbors specified by a range.
+        /** The value_type of an NodeInputIterator should be a size_t, corresponding to
+         *  the indices of nodes of type 2 that should become neighbors of the added node.
+         *  For improved efficiency, the size of the range may be specified by sizeHint.
+         */
+        template <typename NodeInputIterator>
+        void add1( NodeInputIterator begin, NodeInputIterator end, size_t sizeHint = 0 ) {
+            Neighbors nbs1new;
+            nbs1new.reserve( sizeHint );
+            size_t iter = 0;
+            for( NodeInputIterator it = begin; it != end; ++it ) {
+                assert( *it < nr2() );
+                Neighbor nb1new( iter, *it, nb2(*it).size() );
+                Neighbor nb2new( nb2(*it).size(), nr1(), iter++ );
+                nbs1new.push_back( nb1new );
+                nb2( *it ).push_back( nb2new );
             }
+            _nb1.push_back( nbs1new );
+        }
 
-            // Calculate _E12ind
-            _E12ind = std::vector<boost::numeric::ublas::mapped_vector<size_t> >(_V1.size(), boost::numeric::ublas::mapped_vector<size_t>(_V2.size()) );            
-            for( size_t k = 0; k < _E12.size(); k++ )
-                _E12ind[_E12[k].first](_E12[k].second) = k;
-                    
+        /// Add node of type 2 with neighbors specified by a range.
+        /** The value_type of an NodeInputIterator should be a size_t, corresponding to
+         *  the indices of nodes of type 1 that should become neighbors of the added node.
+         *  For improved efficiency, the size of the range may be specified by sizeHint.
+         */
+        template <typename NodeInputIterator>
+        void add2( NodeInputIterator begin, NodeInputIterator end, size_t sizeHint = 0 ) {
+            Neighbors nbs2new;
+            nbs2new.reserve( sizeHint );
+            size_t iter = 0;
+            for( NodeInputIterator it = begin; it != end; ++it ) {
+                assert( *it < nr1() );
+                Neighbor nb2new( iter, *it, nb1(*it).size() );
+                Neighbor nb1new( nb1(*it).size(), nr2(), iter++ );
+                nbs2new.push_back( nb2new );
+                nb1( *it ).push_back( nb1new );
+            }
+            _nb2.push_back( nbs2new );
         }
+
+        /// Remove node of type 1 and all incident edges.
+        void erase1( size_t n1 );
+
+        /// Remove node of type 2 and all incident edges.
+        void erase2( size_t n2 );
+
+        /// Calculate second-order neighbors (i.e., neighbors of neighbors) of node n1 of type 1.
+        /** If include == true, include n1 itself, otherwise exclude n1.
+         */
+        std::vector<size_t> delta1( size_t n1, bool include = false ) const;
+
+        /// Calculate second-order neighbors (i.e., neighbors of neighbors) of node n2 of type 2.
+        /** If include == true, include n2 itself, otherwise exclude n2.
+         */
+        std::vector<size_t> delta2( size_t n2, bool include = false ) const;
+
+        /// Returns true if the graph is connected
+        bool isConnected() const;
+
+        /// Returns true if the graph is a tree, i.e., if it is singly connected and connected.
+        /** This is equivalent to whether for each pair of vertices in the graph, there exists
+         *  a unique path in the graph that starts at the first and ends at the second vertex.
+         */
+        bool isTree() const;
+
+        /// Stream to output stream os in graphviz .dot syntax
+        void display( std::ostream& os ) const;
+
+        /// Checks internal consistency
+        void check() const;
 };
 
 
+template<typename EdgeInputIterator>
+void BipartiteGraph::create( size_t nr1, size_t nr2, EdgeInputIterator begin, EdgeInputIterator end ) {
+    _nb1.clear();
+    _nb1.resize( nr1 );
+    _nb2.clear();
+    _nb2.resize( nr2 );
+
+    for( EdgeInputIterator e = begin; e != end; e++ ) {
+        // Each edge yields a neighbor pair
+        Neighbor nb_1( _nb1[e->first].size(), e->second, _nb2[e->second].size() );
+        Neighbor nb_2( nb_1.dual, e->first, nb_1.iter );
+        _nb1[e->first].push_back( nb_1 );
+        _nb2[e->second].push_back( nb_2 );
+    }
+}
+
+
 } // end of namespace dai