graph.hh

Go to the documentation of this file.

// SPDX-FileCopyrightText: Copyright © DUNE Project contributors, see file LICENSE.md in module root
// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_AMG_GRAPH_HH
#define DUNE_AMG_GRAPH_HH

#include <cstddef>
#include <algorithm>
#include <vector>
#include <cassert>
#include <limits>
#include <dune/common/typetraits.hh>
#include <dune/common/iteratorfacades.hh>
#include <dune/istl/istlexception.hh>
#include <dune/common/propertymap.hh>

namespace Dune
{
  namespace Amg
  {
    template<class M>
    class MatrixGraph
    {
    public:
      typedef M Matrix;

      typedef typename std::remove_const<M>::type MutableMatrix;

      typedef typename M::block_type Weight;

      typedef typename M::size_type VertexDescriptor;

      typedef std::ptrdiff_t EdgeDescriptor;

      enum {
        /*
         * @brief Whether Matrix is mutable.
         */
        mutableMatrix = std::is_same<M, typename std::remove_const<M>::type>::value
      };


      template<class C>
      class EdgeIteratorT
      {

      public:
        typedef typename std::remove_const<C>::type MutableContainer;
        typedef const typename std::remove_const<C>::type ConstContainer;

        friend class EdgeIteratorT<MutableContainer>;
        friend class EdgeIteratorT<ConstContainer>;

        enum {
          isMutable = std::is_same<C, MutableContainer>::value
        };

        typedef typename std::conditional<isMutable && C::mutableMatrix,typename Matrix::row_type::Iterator,
            typename Matrix::row_type::ConstIterator>::type
        ColIterator;

        typedef typename std::conditional<isMutable && C::mutableMatrix,typename M::block_type,
            const typename M::block_type>::type
        Weight;

        EdgeIteratorT(const VertexDescriptor& source, const ColIterator& block,
                      const ColIterator& end, const EdgeDescriptor& edge);

        EdgeIteratorT(const ColIterator& block);

        template<class C1>
        EdgeIteratorT(const EdgeIteratorT<C1>& other);

        typedef typename std::conditional<std::is_same<C, typename std::remove_const<C>::type>::value && C::mutableMatrix,
            typename M::block_type, const typename M::block_type>::type
        WeightType;

        WeightType& weight() const;

        EdgeIteratorT<C>& operator++();

        bool operator!=(const EdgeIteratorT<typename std::remove_const<C>::type>& other) const;

        bool operator!=(const EdgeIteratorT<const typename std::remove_const<C>::type>& other) const;

        bool operator==(const EdgeIteratorT<typename std::remove_const<C>::type>& other) const;

        bool operator==(const EdgeIteratorT<const typename std::remove_const<C>::type>& other) const;

        VertexDescriptor target() const;

        VertexDescriptor source() const;

        const EdgeDescriptor& operator*() const;

        const EdgeDescriptor* operator->() const;

      private:
        VertexDescriptor source_;
        ColIterator block_;
        /***
         * @brief The column iterator positioned at the end of the row
         * of vertex source_
         */
        ColIterator blockEnd_;
        EdgeDescriptor edge_;
      };

      template<class C>
      class VertexIteratorT
      {
      public:
        typedef typename std::remove_const<C>::type MutableContainer;
        typedef const typename std::remove_const<C>::type ConstContainer;

        friend class VertexIteratorT<MutableContainer>;
        friend class VertexIteratorT<ConstContainer>;

        enum {
          isMutable = std::is_same<C, MutableContainer>::value
        };

        explicit VertexIteratorT(C* graph, const VertexDescriptor& current);

        explicit VertexIteratorT(const VertexDescriptor& current);

        VertexIteratorT(const VertexIteratorT<MutableContainer>& other);

        VertexIteratorT<C>& operator++();

        bool operator!=(const VertexIteratorT<ConstContainer>& other) const;

        bool operator==(const VertexIteratorT<ConstContainer>& other) const;

        bool operator!=(const VertexIteratorT<MutableContainer>& other) const;

        bool operator==(const VertexIteratorT<MutableContainer>& other) const;

        typedef typename std::conditional<std::is_same<C, typename std::remove_const<C>::type>::value && C::mutableMatrix,
            typename M::block_type, const typename M::block_type>::type
        WeightType;
        WeightType& weight() const;

        const VertexDescriptor& operator*() const;

        EdgeIteratorT<C> begin() const;

        EdgeIteratorT<C> end() const;

      private:
        C* graph_;
        VertexDescriptor current_;
      };

      typedef EdgeIteratorT<const MatrixGraph<Matrix> > ConstEdgeIterator;

      typedef EdgeIteratorT<MatrixGraph<Matrix> > EdgeIterator;

      typedef VertexIteratorT<const MatrixGraph<Matrix> > ConstVertexIterator;

      typedef VertexIteratorT<MatrixGraph<Matrix> > VertexIterator;

      MatrixGraph(Matrix& matrix);

      ~MatrixGraph();

      VertexIterator begin();

      VertexIterator end();

      ConstVertexIterator begin() const;

      ConstVertexIterator end() const;

      EdgeIterator beginEdges(const VertexDescriptor& source);

      EdgeIterator endEdges(const VertexDescriptor& source);


      ConstEdgeIterator beginEdges(const VertexDescriptor& source) const;

      ConstEdgeIterator endEdges(const VertexDescriptor& source) const;

      Matrix& matrix();

      const Matrix& matrix() const;

      std::size_t noVertices() const;

      VertexDescriptor maxVertex() const;

      std::size_t noEdges() const;

      EdgeDescriptor findEdge(const VertexDescriptor& source,
                              const VertexDescriptor& target) const;

    private:
      Matrix& matrix_;
      EdgeDescriptor* start_;
      MatrixGraph(const MatrixGraph&);

    };

    template<class G, class T>
    class SubGraph
    {
    public:
      typedef G Graph;

      typedef T Excluded;

      typedef typename Graph::VertexDescriptor VertexDescriptor;

      typedef VertexDescriptor* EdgeDescriptor;

      class EdgeIndexMap
      {
      public:
        typedef ReadablePropertyMapTag Category;

        EdgeIndexMap(const EdgeDescriptor& firstEdge)
          : firstEdge_(firstEdge)
        {}

        EdgeIndexMap(const EdgeIndexMap& emap)
          : firstEdge_(emap.firstEdge_)
        {}

        std::size_t operator[](const EdgeDescriptor& edge) const
        {
          return edge-firstEdge_;
        }
      private:
        EdgeDescriptor firstEdge_;
        EdgeIndexMap()
        {}
      };

      EdgeIndexMap getEdgeIndexMap();

      class EdgeIterator : public RandomAccessIteratorFacade<EdgeIterator,const EdgeDescriptor>
      {
      public:
        explicit EdgeIterator(const VertexDescriptor& source, const EdgeDescriptor& edge);

        explicit EdgeIterator(const EdgeDescriptor& edge);

        bool equals(const EdgeIterator& other) const;

        EdgeIterator& increment();

        EdgeIterator& decrement();

        EdgeIterator& advance(std::ptrdiff_t n);

        const EdgeDescriptor& dereference() const;

        const VertexDescriptor& target() const;

        const VertexDescriptor& source() const;

        std::ptrdiff_t distanceTo(const EdgeIterator& other) const;

      private:
        VertexDescriptor source_;
        EdgeDescriptor edge_;
      };

      class VertexIterator
        : public ForwardIteratorFacade<VertexIterator,const VertexDescriptor>
      {
      public:
        explicit VertexIterator(const SubGraph<G,T>* graph, const VertexDescriptor& current,
                                const VertexDescriptor& end);


        explicit VertexIterator(const VertexDescriptor& current);

        VertexIterator& increment();

        bool equals(const VertexIterator& other) const;

        const VertexDescriptor& dereference() const;

        EdgeIterator begin() const;

        EdgeIterator end() const;

      private:
        const SubGraph<Graph,T>* graph_;
        VertexDescriptor current_;
        VertexDescriptor end_;
      };

      typedef EdgeIterator ConstEdgeIterator;

      typedef VertexIterator ConstVertexIterator;

      ConstVertexIterator begin() const;

      ConstVertexIterator end() const;

      ConstEdgeIterator beginEdges(const VertexDescriptor& source) const;

      ConstEdgeIterator endEdges(const VertexDescriptor& source) const;

      std::size_t noVertices() const;

      VertexDescriptor maxVertex() const;

      std::size_t noEdges() const;
      EdgeDescriptor findEdge(const VertexDescriptor& source,
                                    const VertexDescriptor& target) const;
      SubGraph(const Graph& graph, const T& excluded);

      ~SubGraph();

    private:
      const T& excluded_;
      std::size_t noVertices_;
      VertexDescriptor endVertex_;
      int noEdges_;
      VertexDescriptor maxVertex_;
      VertexDescriptor* edges_;
      std::ptrdiff_t* start_;
      std::ptrdiff_t* end_;
      SubGraph(const SubGraph&)
      {}
    };


    template<class G, class VP, class VM=IdentityMap>
    class VertexPropertiesGraph
    {
    public:
      typedef G Graph;

      typedef typename Graph::VertexDescriptor VertexDescriptor;

      typedef typename Graph::EdgeDescriptor EdgeDescriptor;

      typedef VP VertexProperties;

      typedef VM VertexMap;

      typedef typename Graph::EdgeIterator EdgeIterator;

      typedef typename Graph::ConstEdgeIterator ConstEdgeIterator;

      EdgeIterator beginEdges(const VertexDescriptor& source);

      EdgeIterator endEdges(const VertexDescriptor& source);

      ConstEdgeIterator beginEdges(const VertexDescriptor& source) const;

      ConstEdgeIterator endEdges(const VertexDescriptor& source) const;


      template<class C>
      class VertexIteratorT
        : public std::conditional<std::is_same<typename std::remove_const<C>::type,
                  C>::value,
              typename Graph::VertexIterator,
              typename Graph::ConstVertexIterator>::type
      {
        friend class VertexIteratorT<const typename std::remove_const<C>::type>;
        friend class VertexIteratorT<typename std::remove_const<C>::type>;
      public:
        typedef typename std::conditional<std::is_same<typename std::remove_const<C>::type,
                C>::value,
            typename Graph::VertexIterator,
            typename Graph::ConstVertexIterator>::type
        Father;

        typedef typename std::conditional<std::is_same<typename std::remove_const<C>::type,
                C>::value,
            typename Graph::EdgeIterator,
            typename Graph::ConstEdgeIterator>::type
        EdgeIterator;

        explicit VertexIteratorT(const Father& iter,
                                 C* graph);


        explicit VertexIteratorT(const Father& iter);

        template<class C1>
        VertexIteratorT(const VertexIteratorT<C1>& other);

        typename std::conditional<std::is_same<C,typename std::remove_const<C>::type>::value,
            VertexProperties&,
            const VertexProperties&>::type
        properties() const;

        EdgeIterator begin() const;

        EdgeIterator end() const;

      private:
        C* graph_;
      };

      typedef VertexIteratorT<VertexPropertiesGraph<Graph,
              VertexProperties,VM> > VertexIterator;

      typedef VertexIteratorT<const VertexPropertiesGraph<Graph,
              VertexProperties,VM> > ConstVertexIterator;

      VertexIterator begin();

      VertexIterator end();

      ConstVertexIterator begin() const;

      ConstVertexIterator end() const;

      VertexProperties& getVertexProperties(const VertexDescriptor& vertex);

      const VertexProperties& getVertexProperties(const VertexDescriptor& vertex) const;

      const Graph& graph() const;

      std::size_t noVertices() const;

      std::size_t noEdges() const;

      VertexDescriptor maxVertex() const;

      VertexPropertiesGraph(Graph& graph, const VertexMap vmap=VertexMap());

    private:
      VertexPropertiesGraph(const VertexPropertiesGraph&)
      {}

      Graph& graph_;
      VertexMap vmap_;
      std::vector<VertexProperties> vertexProperties_;

    };

    template<class G, class VP, class EP, class VM=IdentityMap, class EM=IdentityMap>
    class PropertiesGraph
    {
    public:
      typedef G Graph;

      typedef typename Graph::VertexDescriptor VertexDescriptor;

      typedef typename Graph::EdgeDescriptor EdgeDescriptor;

      typedef VP VertexProperties;

      typedef VM VertexMap;

      typedef EP EdgeProperties;


      typedef EM EdgeMap;

      template<class C>
      class EdgeIteratorT
        :  public std::conditional<std::is_same<typename std::remove_const<C>::type,
                  C>::value,
              typename Graph::EdgeIterator,
              typename Graph::ConstEdgeIterator>::type
      {

        friend class EdgeIteratorT<const typename std::remove_const<C>::type>;
        friend class EdgeIteratorT<typename std::remove_const<C>::type>;
      public:
        typedef typename std::conditional<std::is_same<typename std::remove_const<C>::type,
                C>::value,
            typename Graph::EdgeIterator,
            typename Graph::ConstEdgeIterator>::type
        Father;

        explicit EdgeIteratorT(const Father& iter,
                               C* graph);

        explicit EdgeIteratorT(const Father& iter);

        template<class C1>
        EdgeIteratorT(const EdgeIteratorT<C1>& other);

        typename std::conditional<std::is_same<C,typename std::remove_const<C>::type>::value,
            EdgeProperties&,
            const EdgeProperties&>::type
        properties() const;

      private:
        C* graph_;
      };

      typedef EdgeIteratorT<PropertiesGraph<Graph,
              VertexProperties,
              EdgeProperties,VM,EM> > EdgeIterator;

      typedef EdgeIteratorT<const PropertiesGraph<Graph,
              VertexProperties,
              EdgeProperties,VM,EM> > ConstEdgeIterator;

      EdgeIterator beginEdges(const VertexDescriptor& source);

      EdgeIterator endEdges(const VertexDescriptor& source);

      ConstEdgeIterator beginEdges(const VertexDescriptor& source) const;

      ConstEdgeIterator endEdges(const VertexDescriptor& source) const;


      template<class C>
      class VertexIteratorT
        : public std::conditional<std::is_same<typename std::remove_const<C>::type,
                  C>::value,
              typename Graph::VertexIterator,
              typename Graph::ConstVertexIterator>::type
      {
        friend class VertexIteratorT<const typename std::remove_const<C>::type>;
        friend class VertexIteratorT<typename std::remove_const<C>::type>;
      public:
        typedef typename std::conditional<std::is_same<typename std::remove_const<C>::type,
                C>::value,
            typename Graph::VertexIterator,
            typename Graph::ConstVertexIterator>::type
        Father;

        explicit VertexIteratorT(const Father& iter,
                                 C* graph);


        explicit VertexIteratorT(const Father& iter);

        template<class C1>
        VertexIteratorT(const VertexIteratorT<C1>& other);

        typename std::conditional<std::is_same<C,typename std::remove_const<C>::type>::value,
            VertexProperties&,
            const VertexProperties&>::type
        properties() const;

        EdgeIteratorT<C> begin() const;

        EdgeIteratorT<C> end() const;

      private:
        C* graph_;
      };

      typedef VertexIteratorT<PropertiesGraph<Graph,
              VertexProperties,
              EdgeProperties,VM,EM> > VertexIterator;

      typedef VertexIteratorT<const PropertiesGraph<Graph,
              VertexProperties,
              EdgeProperties,VM,EM> > ConstVertexIterator;

      VertexIterator begin();

      VertexIterator end();

      ConstVertexIterator begin() const;

      ConstVertexIterator end() const;

      VertexProperties& getVertexProperties(const VertexDescriptor& vertex);

      const VertexProperties& getVertexProperties(const VertexDescriptor& vertex) const;

      EdgeDescriptor findEdge(const VertexDescriptor& source,
                              const VertexDescriptor& target)
      {
        return graph_.findEdge(source,target);
      }

      EdgeProperties& getEdgeProperties(const EdgeDescriptor& edge);


      const EdgeProperties& getEdgeProperties(const EdgeDescriptor& edge) const;

      EdgeProperties& getEdgeProperties(const VertexDescriptor& source,
                                        const VertexDescriptor& target);

      const EdgeProperties& getEdgeProperties(const VertexDescriptor& source,
                                              const VertexDescriptor& target) const;

      const Graph& graph() const;

      std::size_t noVertices() const;

      std::size_t noEdges() const;

      VertexDescriptor maxVertex() const;

      PropertiesGraph(Graph& graph, const VertexMap& vmap=VertexMap(),
                      const EdgeMap& emap=EdgeMap());

    private:
      PropertiesGraph(const PropertiesGraph&)
      {}

      Graph& graph_;
      VertexMap vmap_;
      std::vector<VertexProperties> vertexProperties_;
      EdgeMap emap_;
      std::vector<EdgeProperties> edgeProperties_;

    };


    template<typename G>
    class GraphVertexPropertiesSelector
    {
    public:
      typedef G Graph;
      typedef typename G::VertexProperties VertexProperties;
      typedef typename G::VertexDescriptor Vertex;

      GraphVertexPropertiesSelector(G& g)
        : graph_(g)
      {}
      GraphVertexPropertiesSelector()
        : graph_(0)
      {}


      VertexProperties& operator[](const Vertex& vertex) const
      {
        return graph_->getVertexProperties(vertex);
      }
    private:
      Graph* graph_;
    };

    template<typename G>
    class GraphEdgePropertiesSelector
    {
    public:
      typedef G Graph;
      typedef typename G::EdgeProperties EdgeProperties;
      typedef typename G::EdgeDescriptor Edge;

      GraphEdgePropertiesSelector(G& g)
        : graph_(g)
      {}
      GraphEdgePropertiesSelector()
        : graph_(0)
      {}

      EdgeProperties& operator[](const Edge& edge) const
      {
        return graph_->getEdgeProperties(edge);
      }
    private:
      Graph* graph_;
    };


    template<class G, class V>
    int visitNeighbours(const G& graph, const typename G::VertexDescriptor& vertex,
                        V& visitor);

#ifndef DOXYGEN

    template<class M>
    MatrixGraph<M>::MatrixGraph(M& matrix)
      : matrix_(matrix)
    {
      if(matrix_.N()!=matrix_.M())
        DUNE_THROW(ISTLError, "Matrix has to have as many columns as rows!");

      start_ = new EdgeDescriptor[matrix_.N()+1];

      typedef typename M::ConstIterator Iterator;
      start_[matrix_.begin().index()] = 0;

      for(Iterator row=matrix_.begin(); row != matrix_.end(); ++row)
        start_[row.index()+1] = start_[row.index()] + row->size();
    }

    template<class M>
    MatrixGraph<M>::~MatrixGraph()
    {
      delete[] start_;
    }

    template<class M>
    inline std::size_t MatrixGraph<M>::noEdges() const
    {
      return start_[matrix_.N()];
    }

    template<class M>
    inline std::size_t MatrixGraph<M>::noVertices() const
    {
      return matrix_.N();
    }

    template<class M>
    inline typename MatrixGraph<M>::VertexDescriptor MatrixGraph<M>::maxVertex() const
    {
      return matrix_.N()-1;
    }

    template<class M>
    typename MatrixGraph<M>::EdgeDescriptor
    MatrixGraph<M>::findEdge(const VertexDescriptor& source,
                             const VertexDescriptor& target) const
    {
      typename M::ConstColIterator found =matrix_[source].find(target);
      if(found == matrix_[source].end())
        return std::numeric_limits<EdgeDescriptor>::max();
      std::size_t offset = found.offset();
      if(target>source)
        offset--;

      assert(offset<noEdges());

      return start_[source]+offset;
    }


    template<class M>
    inline M& MatrixGraph<M>::matrix()
    {
      return matrix_;
    }

    template<class M>
    inline const M& MatrixGraph<M>::matrix() const
    {
      return matrix_;
    }

    template<class M>
    template<class C>
    MatrixGraph<M>::EdgeIteratorT<C>::EdgeIteratorT(const VertexDescriptor& source, const ColIterator& block,
                                                    const ColIterator& end, const EdgeDescriptor& edge)
      : source_(source), block_(block), blockEnd_(end), edge_(edge)
    {
      if(block_!=blockEnd_ && block_.index() == source_) {
        // This is the edge from the diagonal to the diagonal. Skip it.
        ++block_;
      }
    }

    template<class M>
    template<class C>
    MatrixGraph<M>::EdgeIteratorT<C>::EdgeIteratorT(const ColIterator& block)
      : block_(block)
    {}

    template<class M>
    template<class C>
    template<class C1>
    MatrixGraph<M>::EdgeIteratorT<C>::EdgeIteratorT(const EdgeIteratorT<C1>& other)
      : source_(other.source_), block_(other.block_), blockEnd_(other.blockEnd_), edge_(other.edge_)
    {}


    template<class M>
    template<class C>
    inline typename MatrixGraph<M>::template EdgeIteratorT<C>::WeightType&
    MatrixGraph<M>::EdgeIteratorT<C>::weight() const
    {
      return *block_;
    }

    template<class M>
    template<class C>
    inline typename MatrixGraph<M>::template EdgeIteratorT<C>& MatrixGraph<M>::EdgeIteratorT<C>::operator++()
    {
      ++block_;
      ++edge_;

      if(block_!=blockEnd_ && block_.index() == source_) {
        // This is the edge from the diagonal to the diagonal. Skip it.
        ++block_;
      }

      return *this;
    }

    template<class M>
    template<class C>
    inline bool MatrixGraph<M>::EdgeIteratorT<C>::operator!=(const typename MatrixGraph<M>::template EdgeIteratorT<typename std::remove_const<C>::type>& other) const
    {
      return block_!=other.block_;
    }

    template<class M>
    template<class C>
    inline bool MatrixGraph<M>::EdgeIteratorT<C>::operator!=(const typename MatrixGraph<M>::template EdgeIteratorT<const typename std::remove_const<C>::type>& other) const
    {
      return block_!=other.block_;
    }

    template<class M>
    template<class C>
    inline bool MatrixGraph<M>::EdgeIteratorT<C>::operator==(const typename MatrixGraph<M>::template EdgeIteratorT<typename std::remove_const<C>::type>& other) const
    {
      return block_==other.block_;
    }

    template<class M>
    template<class C>
    inline bool MatrixGraph<M>::EdgeIteratorT<C>::operator==(const typename MatrixGraph<M>::template EdgeIteratorT<const typename std::remove_const<C>::type>& other) const
    {
      return block_==other.block_;
    }

    template<class M>
    template<class C>
    inline typename MatrixGraph<M>::VertexDescriptor MatrixGraph<M>::EdgeIteratorT<C>::target() const
    {
      return block_.index();
    }

    template<class M>
    template<class C>
    inline typename MatrixGraph<M>::VertexDescriptor MatrixGraph<M>::EdgeIteratorT<C>::source() const
    {
      return source_;
    }

    template<class M>
    template<class C>
    inline const typename MatrixGraph<M>::EdgeDescriptor& MatrixGraph<M>::EdgeIteratorT<C>::operator*() const
    {
      return edge_;
    }

    template<class M>
    template<class C>
    inline const typename MatrixGraph<M>::EdgeDescriptor* MatrixGraph<M>::EdgeIteratorT<C>::operator->() const
    {
      return &edge_;
    }

    template<class M>
    template<class C>
    MatrixGraph<M>::VertexIteratorT<C>::VertexIteratorT(C* graph,
                                                        const VertexDescriptor& current)
      : graph_(graph), current_(current)
    {}


    template<class M>
    template<class C>
    MatrixGraph<M>::VertexIteratorT<C>::VertexIteratorT(const VertexDescriptor& current)
      : current_(current)
    {}

    template<class M>
    template<class C>
    MatrixGraph<M>::VertexIteratorT<C>::VertexIteratorT(const VertexIteratorT<MutableContainer>& other)
      : graph_(other.graph_), current_(other.current_)
    {}

    template<class M>
    template<class C>
    inline bool MatrixGraph<M>::VertexIteratorT<C>::operator!=(const VertexIteratorT<MutableContainer>& other) const
    {
      return current_ != other.current_;
    }

    template<class M>
    template<class C>
    inline bool MatrixGraph<M>::VertexIteratorT<C>::operator!=(const VertexIteratorT<ConstContainer>& other) const
    {
      return current_ != other.current_;
    }


    template<class M>
    template<class C>
    inline bool MatrixGraph<M>::VertexIteratorT<C>::operator==(const VertexIteratorT<MutableContainer>& other) const
    {
      return current_ == other.current_;
    }

    template<class M>
    template<class C>
    inline bool MatrixGraph<M>::VertexIteratorT<C>::operator==(const VertexIteratorT<ConstContainer>& other) const
    {
      return current_ == other.current_;
    }

    template<class M>
    template<class C>
    inline typename MatrixGraph<M>::template VertexIteratorT<C>& MatrixGraph<M>::VertexIteratorT<C>::operator++()
    {
      ++current_;
      return *this;
    }

    template<class M>
    template<class C>
    inline typename MatrixGraph<M>::template VertexIteratorT<C>::WeightType&
    MatrixGraph<M>::VertexIteratorT<C>::weight() const
    {
      return graph_->matrix()[current_][current_];
    }

    template<class M>
    template<class C>
    inline const typename MatrixGraph<M>::VertexDescriptor&
    MatrixGraph<M>::VertexIteratorT<C>::operator*() const
    {
      return current_;
    }

    template<class M>
    template<class C>
    inline typename MatrixGraph<M>::template EdgeIteratorT<C>
    MatrixGraph<M>::VertexIteratorT<C>::begin() const
    {
      return graph_->beginEdges(current_);
    }

    template<class M>
    template<class C>
    inline typename MatrixGraph<M>::template EdgeIteratorT<C>
    MatrixGraph<M>::VertexIteratorT<C>::end() const
    {
      return graph_->endEdges(current_);
    }

    template<class M>
    inline typename MatrixGraph<M>::template VertexIteratorT<MatrixGraph<M> >
    MatrixGraph<M>::begin()
    {
      return VertexIterator(this,0);
    }

    template<class M>
    inline typename MatrixGraph<M>::template VertexIteratorT<MatrixGraph<M> >
    MatrixGraph<M>::end()
    {
      return VertexIterator(matrix_.N());
    }


    template<class M>
    inline typename MatrixGraph<M>::template VertexIteratorT<const MatrixGraph<M> >
    MatrixGraph<M>::begin() const
    {
      return ConstVertexIterator(this, 0);
    }

    template<class M>
    inline typename MatrixGraph<M>::template VertexIteratorT<const MatrixGraph<M> >
    MatrixGraph<M>::end() const
    {
      return ConstVertexIterator(matrix_.N());
    }

    template<class M>
    inline typename MatrixGraph<M>::template EdgeIteratorT<MatrixGraph<M> >
    MatrixGraph<M>::beginEdges(const VertexDescriptor& source)
    {
      return EdgeIterator(source, matrix_.operator[](source).begin(),
                          matrix_.operator[](source).end(), start_[source]);
    }

    template<class M>
    inline typename MatrixGraph<M>::template EdgeIteratorT<MatrixGraph<M> >
    MatrixGraph<M>::endEdges(const VertexDescriptor& source)
    {
      return EdgeIterator(matrix_.operator[](source).end());
    }


    template<class M>
    inline typename MatrixGraph<M>::template EdgeIteratorT<const MatrixGraph<M> >
    MatrixGraph<M>::beginEdges(const VertexDescriptor& source) const
    {
      return ConstEdgeIterator(source, matrix_.operator[](source).begin(),
                               matrix_.operator[](source).end(), start_[source]);
    }

    template<class M>
    inline typename MatrixGraph<M>::template EdgeIteratorT<const MatrixGraph<M> >
    MatrixGraph<M>::endEdges(const VertexDescriptor& source) const
    {
      return ConstEdgeIterator(matrix_.operator[](source).end());
    }


    template<class G, class T>
    SubGraph<G,T>::EdgeIterator::EdgeIterator(const VertexDescriptor& source,
                                              const EdgeDescriptor& edge)
      : source_(source), edge_(edge)
    {}


    template<class G, class T>
    SubGraph<G,T>::EdgeIterator::EdgeIterator(const EdgeDescriptor& edge)
      : edge_(edge)
    {}

    template<class G, class T>
    typename SubGraph<G,T>::EdgeIndexMap SubGraph<G,T>::getEdgeIndexMap()
    {
      return EdgeIndexMap(edges_);
    }

    template<class G, class T>
    inline bool SubGraph<G,T>::EdgeIterator::equals(const EdgeIterator & other) const
    {
      return other.edge_==edge_;
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::EdgeIterator& SubGraph<G,T>::EdgeIterator::increment()
    {
      ++edge_;
      return *this;
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::EdgeIterator& SubGraph<G,T>::EdgeIterator::decrement()
    {
      --edge_;
      return *this;
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::EdgeIterator& SubGraph<G,T>::EdgeIterator::advance(std::ptrdiff_t n)
    {
      edge_+=n;
      return *this;
    }
    template<class G, class T>
    inline const typename G::VertexDescriptor& SubGraph<G,T>::EdgeIterator::source() const
    {
      return source_;
    }

    template<class G, class T>
    inline const typename G::VertexDescriptor& SubGraph<G,T>::EdgeIterator::target() const
    {
      return *edge_;
    }


    template<class G, class T>
    inline const typename SubGraph<G,T>::EdgeDescriptor& SubGraph<G,T>::EdgeIterator::dereference() const
    {
      return edge_;
    }

    template<class G, class T>
    inline std::ptrdiff_t SubGraph<G,T>::EdgeIterator::distanceTo(const EdgeIterator & other) const
    {
      return other.edge_-edge_;
    }

    template<class G, class T>
    SubGraph<G,T>::VertexIterator::VertexIterator(const SubGraph<G,T>* graph,
                                                  const VertexDescriptor& current,
                                                  const VertexDescriptor& end)
      : graph_(graph), current_(current), end_(end)
    {
      // Skip excluded vertices
      typedef typename T::const_iterator Iterator;

      for(Iterator vertex = graph_->excluded_.begin();
          current_ != end_ && *vertex;
          ++vertex)
        ++current_;
      assert(current_ == end_ || !graph_->excluded_[current_]);
    }

    template<class G, class T>
    SubGraph<G,T>::VertexIterator::VertexIterator(const VertexDescriptor& current)
      : current_(current)
    {}

    template<class G, class T>
    inline typename SubGraph<G,T>::VertexIterator& SubGraph<G,T>::VertexIterator::increment()
    {
      ++current_;
      //Skip excluded vertices
      while(current_ != end_ && graph_->excluded_[current_])
        ++current_;

      assert(current_ == end_ || !graph_->excluded_[current_]);
      return *this;
    }

    template<class G, class T>
    inline bool SubGraph<G,T>::VertexIterator::equals(const VertexIterator & other) const
    {
      return current_==other.current_;
    }

    template<class G, class T>
    inline const typename G::VertexDescriptor& SubGraph<G,T>::VertexIterator::dereference() const
    {
      return current_;
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::EdgeIterator SubGraph<G,T>::VertexIterator::begin() const
    {
      return graph_->beginEdges(current_);
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::EdgeIterator SubGraph<G,T>::VertexIterator::end() const
    {
      return graph_->endEdges(current_);
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::VertexIterator SubGraph<G,T>::begin() const
    {
      return VertexIterator(this, 0, endVertex_);
    }


    template<class G, class T>
    inline typename SubGraph<G,T>::VertexIterator SubGraph<G,T>::end() const
    {
      return VertexIterator(endVertex_);
    }


    template<class G, class T>
    inline typename SubGraph<G,T>::EdgeIterator SubGraph<G,T>::beginEdges(const VertexDescriptor& source) const
    {
      return EdgeIterator(source, edges_+start_[source]);
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::EdgeIterator SubGraph<G,T>::endEdges(const VertexDescriptor& source) const
    {
      return EdgeIterator(edges_+end_[source]);
    }

    template<class G, class T>
    std::size_t SubGraph<G,T>::noVertices() const
    {
      return noVertices_;
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::VertexDescriptor SubGraph<G,T>::maxVertex() const
    {
      return maxVertex_;
    }

    template<class G, class T>
    inline std::size_t SubGraph<G,T>::noEdges() const
    {
      return noEdges_;
    }

    template<class G, class T>
    inline typename SubGraph<G,T>::EdgeDescriptor SubGraph<G,T>::findEdge(const VertexDescriptor& source,
                                                                                const VertexDescriptor& target) const
    {
      const EdgeDescriptor edge = std::lower_bound(edges_+start_[source], edges_+end_[source], target);
      if(edge==edges_+end_[source] || *edge!=target)
        return std::numeric_limits<EdgeDescriptor>::max();

      return edge;
    }

    template<class G, class T>
    SubGraph<G,T>::~SubGraph()
    {
      delete[] edges_;
      delete[] end_;
      delete[] start_;
    }

    template<class G, class T>
    SubGraph<G,T>::SubGraph(const G& graph, const T& excluded)
      : excluded_(excluded), noVertices_(0), endVertex_(0), maxVertex_(graph.maxVertex())
    {
      start_ = new std::ptrdiff_t[graph.noVertices()];
      end_ = new std::ptrdiff_t[graph.noVertices()];
      edges_ = new VertexDescriptor[graph.noEdges()];

      VertexDescriptor* edge=edges_;

      // Cater for the case that there are no vertices.
      // Otherwise endVertex_ will get 1 below.
      if ( graph.noVertices() == 0)
        return;

      typedef typename Graph::ConstVertexIterator Iterator;
      Iterator endVertex=graph.end();

      for(Iterator vertex = graph.begin(); vertex != endVertex; ++vertex)
        if(excluded_[*vertex])
          start_[*vertex]=end_[*vertex]=-1;
        else{
          ++noVertices_;
          endVertex_ = std::max(*vertex, endVertex_);

          start_[*vertex] = edge-edges_;

          auto endEdge = vertex.end();

          for(auto iter=vertex.begin(); iter!= endEdge; ++iter)
            if(!excluded[iter.target()]) {
              *edge = iter.target();
              ++edge;
            }

          end_[*vertex] = edge - edges_;

          // Sort the edges
          std::sort(edges_+start_[*vertex], edge);
        }
      noEdges_ = edge-edges_;
      ++endVertex_;
    }

    template<class G, class V, class VM>
    inline std::size_t VertexPropertiesGraph<G,V,VM>::noEdges() const
    {
      return graph_.noEdges();
    }

    template<class G, class V, class VM>
    inline typename VertexPropertiesGraph<G,V,VM>::EdgeIterator
    VertexPropertiesGraph<G,V,VM>::beginEdges(const VertexDescriptor& source)
    {
      return graph_.beginEdges(source);
    }

    template<class G, class V, class VM>
    inline typename VertexPropertiesGraph<G,V,VM>::EdgeIterator
    VertexPropertiesGraph<G,V,VM>::endEdges(const VertexDescriptor& source)
    {
      return graph_.endEdges(source);
    }

    template<class G, class V, class VM>
    typename VertexPropertiesGraph<G,V,VM>::ConstEdgeIterator
    inline VertexPropertiesGraph<G,V,VM>::beginEdges(const VertexDescriptor& source) const
    {
      return graph_.beginEdges(source);
    }

    template<class G, class V, class VM>
    typename VertexPropertiesGraph<G,V,VM>::ConstEdgeIterator
    VertexPropertiesGraph<G,V,VM>::endEdges(const VertexDescriptor& source) const
    {
      return graph_.endEdges(source);
    }

    template<class G, class V, class VM>
    template<class C>
    VertexPropertiesGraph<G,V,VM>::VertexIteratorT<C>
    ::VertexIteratorT(const Father& iter,
                      C* graph)
      : Father(iter), graph_(graph)
    {}

    template<class G, class V, class VM>
    template<class C>
    VertexPropertiesGraph<G,V,VM>::VertexIteratorT<C>
    ::VertexIteratorT(const Father& iter)
      : Father(iter)
    {}

    template<class G, class V, class VM>
    template<class C>
    template<class C1>
    VertexPropertiesGraph<G,V,VM>::VertexIteratorT<C>
    ::VertexIteratorT(const VertexIteratorT<C1>& other)
      : Father(other), graph_(other.graph_)
    {}

    template<class G, class V, class VM>
    template<class C>
    typename std::conditional<std::is_same<C,typename std::remove_const<C>::type>::value,
        V&, const V&>::type
    inline VertexPropertiesGraph<G,V,VM>::VertexIteratorT<C>::properties() const
    {
      return graph_->getVertexProperties(Father::operator*());
    }

    template<class G, class V, class VM>
    template<class C>
    typename std::conditional<std::is_same<typename std::remove_const<C>::type,
            C>::value,
        typename G::EdgeIterator,
        typename G::ConstEdgeIterator>::type
    inline VertexPropertiesGraph<G,V,VM>::VertexIteratorT<C>::begin() const
    {
      return graph_->beginEdges(Father::operator*());
    }

    template<class G, class V, class VM>
    template<class C>
    typename std::conditional<std::is_same<typename std::remove_const<C>::type,
            C>::value,
        typename G::EdgeIterator,
        typename G::ConstEdgeIterator>::type
    inline VertexPropertiesGraph<G,V,VM>::VertexIteratorT<C>::end() const
    {
      return graph_->endEdges(Father::operator*());
    }

    template<class G, class V, class VM>
    inline typename VertexPropertiesGraph<G,V,VM>::VertexIterator VertexPropertiesGraph<G,V,VM>::begin()
    {
      return VertexIterator(graph_.begin(), this);
    }

    template<class G, class V, class VM>
    inline typename VertexPropertiesGraph<G,V,VM>::VertexIterator VertexPropertiesGraph<G,V,VM>::end()
    {
      return VertexIterator(graph_.end());
    }


    template<class G, class V, class VM>
    inline typename VertexPropertiesGraph<G,V,VM>::ConstVertexIterator VertexPropertiesGraph<G,V,VM>::begin() const
    {
      return ConstVertexIterator(graph_.begin(), this);
    }

    template<class G, class V, class VM>
    inline typename VertexPropertiesGraph<G,V,VM>::ConstVertexIterator VertexPropertiesGraph<G,V,VM>::end() const
    {
      return ConstVertexIterator(graph_.end());
    }

    template<class G, class V, class VM>
    inline V& VertexPropertiesGraph<G,V,VM>::getVertexProperties(const VertexDescriptor& vertex)
    {
      return vertexProperties_[vmap_[vertex]];
    }

    template<class G, class V, class VM>
    inline const V& VertexPropertiesGraph<G,V,VM>::getVertexProperties(const VertexDescriptor& vertex) const
    {
      return vertexProperties_[vmap_[vertex]];
    }

    template<class G, class V, class VM>
    inline const G& VertexPropertiesGraph<G,V,VM>::graph() const
    {
      return graph_;
    }

    template<class G, class V, class VM>
    inline std::size_t VertexPropertiesGraph<G,V,VM>::noVertices() const
    {
      return graph_.noVertices();
    }


    template<class G, class V, class VM>
    inline typename VertexPropertiesGraph<G,V,VM>::VertexDescriptor VertexPropertiesGraph<G,V,VM>::maxVertex() const
    {
      return graph_.maxVertex();
    }

    template<class G, class V, class VM>
    VertexPropertiesGraph<G,V,VM>::VertexPropertiesGraph(Graph& graph, const VM vmap)
      : graph_(graph), vmap_(vmap), vertexProperties_(vmap_[graph_.maxVertex()+1], V())
    {}

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    PropertiesGraph<G,V,E,VM,EM>::EdgeIteratorT<C>::EdgeIteratorT(const Father& iter,
                                                                  C* graph)
      : Father(iter), graph_(graph)
    {}

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    PropertiesGraph<G,V,E,VM,EM>::EdgeIteratorT<C>::EdgeIteratorT(const Father& iter)
      : Father(iter)
    {}

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    template<class C1>
    PropertiesGraph<G,V,E,VM,EM>::EdgeIteratorT<C>::EdgeIteratorT(const EdgeIteratorT<C1>& other)
      : Father(other), graph_(other.graph_)
    {}


    template<class G, class V, class E, class VM, class EM>
    inline std::size_t PropertiesGraph<G,V,E,VM,EM>::noEdges() const
    {
      return graph_.noEdges();
    }

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    inline typename std::conditional<std::is_same<C,typename std::remove_const<C>::type>::value,E&,const E&>::type
    PropertiesGraph<G,V,E,VM,EM>::EdgeIteratorT<C>::properties() const
    {
      return graph_->getEdgeProperties(Father::operator*());
    }

    template<class G, class V, class E, class VM, class EM>
    inline typename PropertiesGraph<G,V,E,VM,EM>::EdgeIterator
    PropertiesGraph<G,V,E,VM,EM>::beginEdges(const VertexDescriptor& source)
    {
      return EdgeIterator(graph_.beginEdges(source), this);
    }

    template<class G, class V, class E, class VM, class EM>
    inline typename PropertiesGraph<G,V,E,VM,EM>::EdgeIterator
    PropertiesGraph<G,V,E,VM,EM>::endEdges(const VertexDescriptor& source)
    {
      return EdgeIterator(graph_.endEdges(source));
    }

    template<class G, class V, class E, class VM, class EM>
    typename PropertiesGraph<G,V,E,VM,EM>::ConstEdgeIterator
    inline PropertiesGraph<G,V,E,VM,EM>::beginEdges(const VertexDescriptor& source) const
    {
      return ConstEdgeIterator(graph_.beginEdges(source), this);
    }

    template<class G, class V, class E, class VM, class EM>
    typename PropertiesGraph<G,V,E,VM,EM>::ConstEdgeIterator
    PropertiesGraph<G,V,E,VM,EM>::endEdges(const VertexDescriptor& source) const
    {
      return ConstEdgeIterator(graph_.endEdges(source));
    }

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    PropertiesGraph<G,V,E,VM,EM>::VertexIteratorT<C>
    ::VertexIteratorT(const Father& iter,
                      C* graph)
      : Father(iter), graph_(graph)
    {}

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    PropertiesGraph<G,V,E,VM,EM>::VertexIteratorT<C>
    ::VertexIteratorT(const Father& iter)
      : Father(iter)
    {}

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    template<class C1>
    PropertiesGraph<G,V,E,VM,EM>::VertexIteratorT<C>
    ::VertexIteratorT(const VertexIteratorT<C1>& other)
      : Father(other), graph_(other.graph_)
    {}

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    inline typename std::conditional<std::is_same<C,typename std::remove_const<C>::type>::value,
        V&, const V&>::type
    PropertiesGraph<G,V,E,VM,EM>::VertexIteratorT<C>::properties() const
    {
      return graph_->getVertexProperties(Father::operator*());
    }

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    inline typename PropertiesGraph<G,V,E,VM,EM>::template EdgeIteratorT<C>
    PropertiesGraph<G,V,E,VM,EM>::VertexIteratorT<C>::begin() const
    {
      return graph_->beginEdges(Father::operator*());
    }

    template<class G, class V, class E, class VM, class EM>
    template<class C>
    inline typename PropertiesGraph<G,V,E,VM,EM>::template EdgeIteratorT<C>
    PropertiesGraph<G,V,E,VM,EM>::VertexIteratorT<C>::end() const
    {
      return graph_->endEdges(Father::operator*());
    }

    template<class G, class V, class E, class VM, class EM>
    inline typename PropertiesGraph<G,V,E,VM,EM>::VertexIterator PropertiesGraph<G,V,E,VM,EM>::begin()
    {
      return VertexIterator(graph_.begin(), this);
    }

    template<class G, class V, class E, class VM, class EM>
    inline typename PropertiesGraph<G,V,E,VM,EM>::VertexIterator PropertiesGraph<G,V,E,VM,EM>::end()
    {
      return VertexIterator(graph_.end());
    }


    template<class G, class V, class E, class VM, class EM>
    inline typename PropertiesGraph<G,V,E,VM,EM>::ConstVertexIterator PropertiesGraph<G,V,E,VM,EM>::begin() const
    {
      return ConstVertexIterator(graph_.begin(), this);
    }

    template<class G, class V, class E, class VM, class EM>
    inline typename PropertiesGraph<G,V,E,VM,EM>::ConstVertexIterator PropertiesGraph<G,V,E,VM,EM>::end() const
    {
      return ConstVertexIterator(graph_.end());
    }

    template<class G, class V, class E, class VM, class EM>
    inline V& PropertiesGraph<G,V,E,VM,EM>::getVertexProperties(const VertexDescriptor& vertex)
    {
      return vertexProperties_[vmap_[vertex]];
    }

    template<class G, class V, class E, class VM, class EM>
    inline const V& PropertiesGraph<G,V,E,VM,EM>::getVertexProperties(const VertexDescriptor& vertex) const
    {
      return vertexProperties_[vmap_[vertex]];
    }

    template<class G, class V, class E, class VM, class EM>
    inline E& PropertiesGraph<G,V,E,VM,EM>::getEdgeProperties(const EdgeDescriptor& edge)
    {
      return edgeProperties_[emap_[edge]];
    }

    template<class G, class V, class E, class VM, class EM>
    inline const E& PropertiesGraph<G,V,E,VM,EM>::getEdgeProperties(const EdgeDescriptor& edge) const
    {
      return edgeProperties_[emap_[edge]];
    }

    template<class G, class V, class E, class VM, class EM>
    inline E& PropertiesGraph<G,V,E,VM,EM>::getEdgeProperties(const VertexDescriptor& source,
                                                              const VertexDescriptor& target)
    {
      return getEdgeProperties(graph_.findEdge(source,target));
    }

    template<class G, class V, class E, class VM, class EM>
    inline const E& PropertiesGraph<G,V,E,VM,EM>::getEdgeProperties(const VertexDescriptor& source,
                                                                    const VertexDescriptor& target) const
    {
      return getEdgeProperties(graph_.findEdge(source,target));
    }

    template<class G, class V, class E, class VM, class EM>
    inline const G& PropertiesGraph<G,V,E,VM,EM>::graph() const
    {
      return graph_;
    }

    template<class G, class V, class E, class VM, class EM>
    inline std::size_t PropertiesGraph<G,V,E,VM,EM>::noVertices() const
    {
      return graph_.noVertices();
    }


    template<class G, class V, class E, class VM, class EM>
    inline typename PropertiesGraph<G,V,E,VM,EM>::VertexDescriptor PropertiesGraph<G,V,E,VM,EM>::maxVertex() const
    {
      return graph_.maxVertex();
    }

    template<class G, class V, class E, class VM, class EM>
    PropertiesGraph<G,V,E,VM,EM>::PropertiesGraph(Graph& graph, const VM& vmap, const EM& emap)
      : graph_(graph), vmap_(vmap), vertexProperties_(vmap_[graph_.maxVertex()+1], V()),
        emap_(emap), edgeProperties_(graph_.noEdges(), E())
    {}

    template<class G, class V>
    inline int visitNeighbours(const G& graph, const typename G::VertexDescriptor& vertex,
                               V& visitor)
    {
      typedef typename G::ConstEdgeIterator iterator;
      const iterator end = graph.endEdges(vertex);
      int noNeighbours=0;
      for(iterator edge = graph.beginEdges(vertex); edge != end; ++edge, ++noNeighbours)
        visitor(edge);
      return noNeighbours;
    }

#endif // DOXYGEN

  }
}
#endif