Index: lemon/vf2.h
===================================================================
--- lemon/vf2.h	(revision 1370)
+++ lemon/vf2.h	(revision 1405)
@@ -3,5 +3,5 @@
  * This file is a part of LEMON, a generic C++ optimization library.
  *
- * Copyright (C) 2015
+ * Copyright (C) 2015-2017
  * EMAXA Kutato-fejleszto Kft. (EMAXA Research Ltd.)
  *
@@ -27,20 +27,16 @@
 #include <lemon/dfs.h>
 #include <lemon/bfs.h>
+#include <lemon/bits/vf2_internals.h>
 
 #include <vector>
-#include <set>
-
-namespace lemon
-{
-  namespace bits
-  {
-    namespace vf2
-    {
-      class AlwaysEq
-      {
+
+namespace lemon {
+  namespace bits {
+    namespace vf2 {
+
+      class AlwaysEq {
       public:
         template<class T1, class T2>
-        bool operator()(T1, T2) const
-        {
+        bool operator()(T1&, T2&) const {
           return true;
         }
@@ -48,19 +44,18 @@
 
       template<class M1, class M2>
-      class MapEq
-      {
+      class MapEq{
         const M1 &_m1;
         const M2 &_m2;
       public:
-        MapEq(const M1 &m1, const M2 &m2) : _m1(m1), _m2(m2) {}
-        bool operator()(typename M1::Key k1, typename M2::Key k2) const
-        {
+        MapEq(const M1 &m1, const M2 &m2) : _m1(m1), _m2(m2) { }
+        bool operator()(typename M1::Key k1, typename M2::Key k2) const {
           return _m1[k1] == _m2[k2];
         }
       };
 
+
+
       template <class G>
-      class DfsLeaveOrder : public DfsVisitor<G>
-      {
+      class DfsLeaveOrder : public DfsVisitor<G> {
         const G &_g;
         std::vector<typename G::Node> &_order;
@@ -68,8 +63,6 @@
       public:
         DfsLeaveOrder(const G &g, std::vector<typename G::Node> &order)
-          : i(countNodes(g)), _g(g), _order(order)
-        {}
-        void leave(const typename G::Node &node)
-        {
+          : i(countNodes(g)), _g(g), _order(order) { }
+        void leave(const typename G::Node &node) {
           _order[--i]=node;
         }
@@ -77,6 +70,5 @@
 
       template <class G>
-      class BfsLeaveOrder : public BfsVisitor<G>
-      {
+      class BfsLeaveOrder : public BfsVisitor<G> {
         int i;
         const G &_g;
@@ -84,8 +76,7 @@
       public:
         BfsLeaveOrder(const G &g, std::vector<typename G::Node> &order)
-          : i(0), _g(g), _order(order)
-        {}
-        void process(const typename G::Node &node)
-        {
+          : i(0), _g(g), _order(order){
+        }
+        void process(const typename G::Node &node) {
           _order[i++]=node;
         }
@@ -94,26 +85,4 @@
   }
 
-  ///Graph mapping types.
-
-  ///\ingroup graph_isomorphism
-  ///The \ref Vf2 "VF2" algorithm is capable of finding different kind of
-  ///embeddings, this enum specifies its type.
-  ///
-  ///See \ref graph_isomorphism for a more detailed description.
-  enum Vf2MappingType {
-    /// Subgraph isomorphism
-    SUBGRAPH = 0,
-    /// Induced subgraph isomorphism
-    INDUCED = 1,
-    /// Graph isomorphism
-
-    /// If the two graph has the same number of nodes, than it is
-    /// equivalent to \ref INDUCED, and if they also have the same
-    /// number of edges, then it is also equivalent to \ref SUBGRAPH.
-    ///
-    /// However, using this setting is faster than the other two
-    /// options.
-    ISOMORPH = 2
-  };
 
   ///%VF2 algorithm class.
@@ -145,10 +114,9 @@
            class NEQ = bits::vf2::AlwaysEq >
 #endif
-  class Vf2
-  {
+  class Vf2 {
     //Current depth in the DFS tree.
     int _depth;
     //Functor with bool operator()(G1::Node,G2::Node), which returns 1
-    //if and only if the 2 nodes are equivalent.
+    //ifff the two nodes are equivalent.
     NEQ _nEq;
 
@@ -157,5 +125,5 @@
     //a node of g2.
     M &_mapping;
-    //order[i] is the node of g1, for which we find a pair if depth=i
+    //order[i] is the node of g1, for which we search a pair if depth=i
     std::vector<typename G1::Node> order;
     //currEdgeIts[i] is an edge iterator, witch is last used in the ith
@@ -164,109 +132,102 @@
     //The small graph.
     const G1 &_g1;
-    //The big graph.
+    //The large graph.
     const G2 &_g2;
-    //lookup tables for cut the searchtree
+    //lookup tables for cutting the searchtree
     typename G1::template NodeMap<int> rNew1t,rInOut1t;
 
-    Vf2MappingType _mapping_type;
+    MappingType _mapping_type;
+
+    bool _deallocMappingAfterUse;
 
     //cut the search tree
-    template<Vf2MappingType MT>
-    bool cut(const typename G1::Node n1,const typename G2::Node n2) const
-    {
+    template<MappingType MT>
+    bool cut(const typename G1::Node n1,const typename G2::Node n2) const {
       int rNew2=0,rInOut2=0;
-      for(typename G2::IncEdgeIt e2(_g2,n2); e2!=INVALID; ++e2)
-        {
-          const typename G2::Node currNode=_g2.oppositeNode(n2,e2);
-          if(_conn[currNode]>0)
-            ++rInOut2;
-          else if(MT!=SUBGRAPH&&_conn[currNode]==0)
-            ++rNew2;
-        }
-      switch(MT)
-        {
-        case INDUCED:
-          return rInOut1t[n1]<=rInOut2&&rNew1t[n1]<=rNew2;
-        case SUBGRAPH:
-          return rInOut1t[n1]<=rInOut2;
-        case ISOMORPH:
-          return rInOut1t[n1]==rInOut2&&rNew1t[n1]==rNew2;
-        default:
-          return false;
-        }
-    }
-
-    template<Vf2MappingType MT>
-    bool feas(const typename G1::Node n1,const typename G2::Node n2)
-    {
+      for(typename G2::IncEdgeIt e2(_g2,n2); e2!=INVALID; ++e2) {
+        const typename G2::Node currNode=_g2.oppositeNode(n2,e2);
+        if(_conn[currNode]>0)
+          ++rInOut2;
+        else if(MT!=SUBGRAPH&&_conn[currNode]==0)
+          ++rNew2;
+      }
+      switch(MT) {
+      case INDUCED:
+        return rInOut1t[n1]<=rInOut2&&rNew1t[n1]<=rNew2;
+      case SUBGRAPH:
+        return rInOut1t[n1]<=rInOut2;
+      case ISOMORPH:
+        return rInOut1t[n1]==rInOut2&&rNew1t[n1]==rNew2;
+      default:
+        return false;
+      }
+    }
+
+    template<MappingType MT>
+    bool feas(const typename G1::Node n1,const typename G2::Node n2) {
       if(!_nEq(n1,n2))
         return 0;
 
-      for(typename G1::IncEdgeIt e1(_g1,n1); e1!=INVALID; ++e1)
-        {
-          const typename G1::Node currNode=_g1.oppositeNode(n1,e1);
-          if(_mapping[currNode]!=INVALID)
-            --_conn[_mapping[currNode]];
-        }
+      for(typename G1::IncEdgeIt e1(_g1,n1); e1!=INVALID; ++e1) {
+        const typename G1::Node& currNode=_g1.oppositeNode(n1,e1);
+        if(_mapping[currNode]!=INVALID)
+          --_conn[_mapping[currNode]];
+      }
       bool isIso=1;
-      for(typename G2::IncEdgeIt e2(_g2,n2); e2!=INVALID; ++e2)
-        {
-          const typename G2::Node currNode=_g2.oppositeNode(n2,e2);
-          if(_conn[currNode]<-1)
-            ++_conn[currNode];
-          else if(MT!=SUBGRAPH&&_conn[currNode]==-1)
-            {
+      for(typename G2::IncEdgeIt e2(_g2,n2); e2!=INVALID; ++e2) {
+        int& connCurrNode = _conn[_g2.oppositeNode(n2,e2)];
+        if(connCurrNode<-1)
+          ++connCurrNode;
+        else if(MT!=SUBGRAPH&&connCurrNode==-1) {
+          isIso=0;
+          break;
+        }
+      }
+
+      for(typename G1::IncEdgeIt e1(_g1,n1); e1!=INVALID; ++e1) {
+        const typename G2::Node& currNodePair=_mapping[_g1.oppositeNode(n1,e1)];
+        int& connCurrNodePair=_conn[currNodePair];
+        if(currNodePair!=INVALID&&connCurrNodePair!=-1) {
+          switch(MT) {
+          case INDUCED:
+          case ISOMORPH:
+            isIso=0;
+            break;
+          case SUBGRAPH:
+            if(connCurrNodePair<-1)
               isIso=0;
-              break;
-            }
-        }
-
-      for(typename G1::IncEdgeIt e1(_g1,n1); e1!=INVALID; ++e1)
-        {
-          const typename G1::Node currNode=_g1.oppositeNode(n1,e1);
-          if(_mapping[currNode]!=INVALID&&_conn[_mapping[currNode]]!=-1)
-            {
-              switch(MT)
-                {
-                case INDUCED:
-                case ISOMORPH:
-                  isIso=0;
-                  break;
-                case SUBGRAPH:
-                  if(_conn[_mapping[currNode]]<-1)
-                    isIso=0;
-                  break;
-                }
-              _conn[_mapping[currNode]]=-1;
-            }
-        }
+            break;
+          }
+          connCurrNodePair=-1;
+        }
+      }
       return isIso&&cut<MT>(n1,n2);
     }
 
-    void addPair(const typename G1::Node n1,const typename G2::Node n2)
-    {
+    void addPair(const typename G1::Node n1,const typename G2::Node n2) {
       _conn[n2]=-1;
       _mapping.set(n1,n2);
-      for(typename G2::IncEdgeIt e2(_g2,n2); e2!=INVALID; ++e2)
-        if(_conn[_g2.oppositeNode(n2,e2)]!=-1)
-          ++_conn[_g2.oppositeNode(n2,e2)];
-    }
-
-    void subPair(const typename G1::Node n1,const typename G2::Node n2)
-    {
+      for(typename G2::IncEdgeIt e2(_g2,n2); e2!=INVALID; ++e2) {
+        int& currConn = _conn[_g2.oppositeNode(n2,e2)];
+        if(currConn!=-1)
+          ++currConn;
+      }
+    }
+
+    void subPair(const typename G1::Node n1,const typename G2::Node n2) {
       _conn[n2]=0;
       _mapping.set(n1,INVALID);
-      for(typename G2::IncEdgeIt e2(_g2,n2); e2!=INVALID; ++e2)
-        {
-          const typename G2::Node currNode=_g2.oppositeNode(n2,e2);
-          if(_conn[currNode]>0)
-            --_conn[currNode];
-          else if(_conn[currNode]==-1)
-            ++_conn[n2];
-        }
-    }
-
-    void setOrder()//we will find pairs for the nodes of g1 in this order
-    {
+      for(typename G2::IncEdgeIt e2(_g2,n2); e2!=INVALID; ++e2) {
+        int& currConn = _conn[_g2.oppositeNode(n2,e2)];
+        if(currConn>0)
+          --currConn;
+        else if(currConn==-1)
+          ++_conn[n2];
+      }
+    }
+
+    void setOrder() {
+      // we will find pairs for the nodes of g1 in this order
+
       // bits::vf2::DfsLeaveOrder<G1> v(_g1,order);
       //   DfsVisit<G1,bits::vf2::DfsLeaveOrder<G1> >dfs(_g1, v);
@@ -279,115 +240,101 @@
     }
 
-    template<Vf2MappingType MT>
-    bool extMatch()
-    {
-      while(_depth>=0)
-        {
-          //there are not nodes in g1, which has not pair in g2.
-          if(_depth==static_cast<int>(order.size()))
-            {
+    template<MappingType MT>
+    bool extMatch() {
+      while(_depth>=0) {
+        //there are not nodes in g1, which has not pair in g2.
+        if(_depth==static_cast<int>(order.size())) {
+          --_depth;
+          return true;
+        }
+        typename G1::Node& nodeOfDepth = order[_depth];
+        const typename G2::Node& pairOfNodeOfDepth = _mapping[nodeOfDepth];
+        typename G2::IncEdgeIt &edgeItOfDepth = currEdgeIts[_depth];
+        //the node of g2, which neighbours are the candidates for
+        //the pair of nodeOfDepth
+        typename G2::Node currPNode;
+        if(edgeItOfDepth==INVALID) {
+          typename G1::IncEdgeIt fstMatchedE(_g1,nodeOfDepth);
+          //if pairOfNodeOfDepth!=INVALID, we dont use
+          //fstMatchedE
+          if(pairOfNodeOfDepth==INVALID)
+            for(; fstMatchedE!=INVALID &&
+                  _mapping[_g1.oppositeNode(nodeOfDepth,
+                                            fstMatchedE)]==INVALID;
+                ++fstMatchedE) ; //find fstMatchedE
+          if(fstMatchedE==INVALID||pairOfNodeOfDepth!=INVALID) {
+            //We found no covered neighbours, this means
+            //the graph is not connected(or _depth==0).  Each
+            //uncovered(and there are some other properties due
+            //to the spec. problem types) node of g2 is
+            //candidate.  We can read the iterator of the last
+            //tried node from the match if it is not the first
+            //try(match[nodeOfDepth]!=INVALID)
+            typename G2::NodeIt n2(_g2);
+            //if it's not the first try
+            if(pairOfNodeOfDepth!=INVALID) {
+              n2=++typename G2::NodeIt(_g2,pairOfNodeOfDepth);
+              subPair(nodeOfDepth,pairOfNodeOfDepth);
+            }
+            for(; n2!=INVALID; ++n2)
+              if(MT!=SUBGRAPH) {
+                if(_conn[n2]==0&&feas<MT>(nodeOfDepth,n2))
+                  break;
+              }
+              else if(_conn[n2]>=0&&feas<MT>(nodeOfDepth,n2))
+                break;
+            // n2 is the next candidate
+            if(n2!=INVALID){
+              addPair(nodeOfDepth,n2);
+              ++_depth;
+            }
+            else // there are no more candidates
               --_depth;
-              return true;
-            }
-          //the node of g2, which neighbours are the candidates for
-          //the pair of order[_depth]
-          typename G2::Node currPNode;
-          if(currEdgeIts[_depth]==INVALID)
-            {
-              typename G1::IncEdgeIt fstMatchedE(_g1,order[_depth]);
-              //if _mapping[order[_depth]]!=INVALID, we dont use
-              //fstMatchedE
-              if(_mapping[order[_depth]]==INVALID)
-                for(; fstMatchedE!=INVALID &&
-                      _mapping[_g1.oppositeNode(order[_depth],
-                                              fstMatchedE)]==INVALID;
-                    ++fstMatchedE) ; //find fstMatchedE
-              if(fstMatchedE==INVALID||_mapping[order[_depth]]!=INVALID)
-                {
-                  //We did not find an covered neighbour, this means
-                  //the graph is not connected(or _depth==0).  Every
-                  //uncovered(and there are some other properties due
-                  //to the spec. problem types) node of g2 is
-                  //candidate.  We can read the iterator of the last
-                  //tryed node from the match if it is not the first
-                  //try(match[order[_depth]]!=INVALID)
-                  typename G2::NodeIt n2(_g2);
-                  //if its not the first try
-                  if(_mapping[order[_depth]]!=INVALID)
-                    {
-                      n2=++typename G2::NodeIt(_g2,_mapping[order[_depth]]);
-                      subPair(order[_depth],_mapping[order[_depth]]);
-                    }
-                  for(; n2!=INVALID; ++n2)
-                    if(MT!=SUBGRAPH&&_conn[n2]==0)
-                      {
-                        if(feas<MT>(order[_depth],n2))
-                          break;
-                      }
-                    else if(MT==SUBGRAPH&&_conn[n2]>=0)
-                      if(feas<MT>(order[_depth],n2))
-                        break;
-                  // n2 is the next candidate
-                  if(n2!=INVALID)
-                    {
-                      addPair(order[_depth],n2);
-                      ++_depth;
-                    }
-                  else // there is no more candidate
-                    --_depth;
-                  continue;
-                }
-              else
-                {
-                  currPNode=_mapping[_g1.oppositeNode(order[_depth],
-                                                      fstMatchedE)];
-                  currEdgeIts[_depth]=typename G2::IncEdgeIt(_g2,currPNode);
-                }
-            }
-          else
-            {
-              currPNode=_g2.oppositeNode(_mapping[order[_depth]],
-                                         currEdgeIts[_depth]);
-              subPair(order[_depth],_mapping[order[_depth]]);
-              ++currEdgeIts[_depth];
-            }
-          for(; currEdgeIts[_depth]!=INVALID; ++currEdgeIts[_depth])
-            {
-              const typename G2::Node currNode =
-                _g2.oppositeNode(currPNode, currEdgeIts[_depth]);
-              //if currNode is uncovered
-              if(_conn[currNode]>0&&feas<MT>(order[_depth],currNode))
-                {
-                  addPair(order[_depth],currNode);
-                  break;
-                }
-            }
-          currEdgeIts[_depth]==INVALID?--_depth:++_depth;
-        }
+            continue;
+          }
+          else {
+            currPNode=_mapping[_g1.oppositeNode(nodeOfDepth,
+                                                fstMatchedE)];
+            edgeItOfDepth=typename G2::IncEdgeIt(_g2,currPNode);
+          }
+        }
+        else {
+          currPNode=_g2.oppositeNode(pairOfNodeOfDepth,
+                                     edgeItOfDepth);
+          subPair(nodeOfDepth,pairOfNodeOfDepth);
+          ++edgeItOfDepth;
+        }
+        for(; edgeItOfDepth!=INVALID; ++edgeItOfDepth) {
+          const typename G2::Node currNode =
+            _g2.oppositeNode(currPNode, edgeItOfDepth);
+          if(_conn[currNode]>0&&feas<MT>(nodeOfDepth,currNode)) {
+            addPair(nodeOfDepth,currNode);
+            break;
+          }
+        }
+        edgeItOfDepth==INVALID?--_depth:++_depth;
+      }
       return false;
     }
 
     //calc. the lookup table for cut the searchtree
-    void setRNew1tRInOut1t()
-    {
+    void setRNew1tRInOut1t() {
       typename G1::template NodeMap<int> tmp(_g1,0);
-      for(unsigned int i=0; i<order.size(); ++i)
-        {
-          tmp[order[i]]=-1;
-          for(typename G1::IncEdgeIt e1(_g1,order[i]); e1!=INVALID; ++e1)
-            {
-              const typename G1::Node currNode=_g1.oppositeNode(order[i],e1);
-              if(tmp[currNode]>0)
-                ++rInOut1t[order[i]];
-              else if(tmp[currNode]==0)
-                ++rNew1t[order[i]];
-            }
-          for(typename G1::IncEdgeIt e1(_g1,order[i]); e1!=INVALID; ++e1)
-            {
-              const typename G1::Node currNode=_g1.oppositeNode(order[i],e1);
-              if(tmp[currNode]!=-1)
-                ++tmp[currNode];
-            }
-        }
+      for(unsigned int i=0; i<order.size(); ++i) {
+        const typename G1::Node& orderI = order[i];
+        tmp[orderI]=-1;
+        for(typename G1::IncEdgeIt e1(_g1,orderI); e1!=INVALID; ++e1) {
+          const typename G1::Node currNode=_g1.oppositeNode(orderI,e1);
+          if(tmp[currNode]>0)
+            ++rInOut1t[orderI];
+          else if(tmp[currNode]==0)
+            ++rNew1t[orderI];
+        }
+        for(typename G1::IncEdgeIt e1(_g1,orderI); e1!=INVALID; ++e1) {
+          const typename G1::Node currNode=_g1.oppositeNode(orderI,e1);
+          if(tmp[currNode]!=-1)
+            ++tmp[currNode];
+        }
+      }
     }
   public:
@@ -400,11 +347,10 @@
     ///\param m \ref concepts::ReadWriteMap "read-write" NodeMap
     ///storing the found mapping.
-    ///\param neq A bool-valued binary functor determinining whether a node is
+    ///\param neq A bool-valued binary functor determining whether a node is
     ///mappable to another. By default it is an always true operator.
-    Vf2(const G1 &g1, const G2 &g2,M &m, const NEQ &neq = NEQ() ) :
+    Vf2(const G1 &g1, const G2 &g2, M &m, const NEQ &neq = NEQ() ) :
       _nEq(neq), _conn(g2,0), _mapping(m), order(countNodes(g1)),
       currEdgeIts(countNodes(g1),INVALID), _g1(g1), _g2(g2), rNew1t(g1,0),
-      rInOut1t(g1,0), _mapping_type(SUBGRAPH)
-    {
+      rInOut1t(g1,0), _mapping_type(SUBGRAPH), _deallocMappingAfterUse(0) {
       _depth=0;
       setOrder();
@@ -414,9 +360,21 @@
     }
 
+    ///Destructor
+
+    ///Destructor.
+    ///
+
+    ~Vf2(){
+      if(_deallocMappingAfterUse)
+        delete &_mapping;
+    }
+
     ///Returns the current mapping type
 
     ///Returns the current mapping type
     ///
-    Vf2MappingType mappingType() const { return _mapping_type; }
+    MappingType mappingType() const {
+      return _mapping_type;
+    }
     ///Sets mapping type
 
@@ -425,11 +383,13 @@
     ///The mapping type is set to \ref SUBGRAPH by default.
     ///
-    ///\sa See \ref Vf2MappingType for the possible values.
-    void mappingType(Vf2MappingType m_type) { _mapping_type = m_type; }
+    ///\sa See \ref MappingType for the possible values.
+    void mappingType(MappingType m_type) {
+      _mapping_type = m_type;
+    }
 
     ///Finds a mapping
 
-    ///It finds a mapping between from g1 into g2 according to the mapping
-    ///type set by \ref mappingType(Vf2MappingType) "mappingType()".
+    ///It finds a mapping from g1 into g2 according to the mapping
+    ///type set by \ref mappingType(MappingType) "mappingType()".
     ///
     ///By subsequent calls, it returns all possible mappings one-by-one.
@@ -437,23 +397,20 @@
     ///\retval true if a mapping is found.
     ///\retval false if there is no (more) mapping.
-    bool find()
-    {
-      switch(_mapping_type)
-        {
-        case SUBGRAPH:
-          return extMatch<SUBGRAPH>();
-        case INDUCED:
-          return extMatch<INDUCED>();
-        case ISOMORPH:
-          return extMatch<ISOMORPH>();
-        default:
-          return false;
-        }
+    bool find() {
+      switch(_mapping_type) {
+      case SUBGRAPH:
+        return extMatch<SUBGRAPH>();
+      case INDUCED:
+        return extMatch<INDUCED>();
+      case ISOMORPH:
+        return extMatch<ISOMORPH>();
+      default:
+        return false;
+      }
     }
   };
 
   template<class G1, class G2>
-  class Vf2WizardBase
-  {
+  class Vf2WizardBase {
   protected:
     typedef G1 Graph1;
@@ -463,13 +420,15 @@
     const G2 &_g2;
 
-    Vf2MappingType _mapping_type;
+    MappingType _mapping_type;
 
     typedef typename G1::template NodeMap<typename G2::Node> Mapping;
     bool _local_mapping;
     void *_mapping;
-    void createMapping()
-    {
+    void createMapping() {
       _mapping = new Mapping(_g1);
     }
+
+    void *myVf2; //used in Vf2Wizard::find
+
 
     typedef bits::vf2::AlwaysEq NodeEq;
@@ -477,6 +436,7 @@
 
     Vf2WizardBase(const G1 &g1,const G2 &g2)
-      : _g1(g1), _g2(g2), _mapping_type(SUBGRAPH), _local_mapping(true) {}
+      : _g1(g1), _g2(g2), _mapping_type(SUBGRAPH), _local_mapping(true) { }
   };
+
 
   /// Auxiliary class for the function-type interface of %VF2 algorithm.
@@ -490,6 +450,5 @@
   /// algorithm.
   template<class TR>
-  class Vf2Wizard : public TR
-  {
+  class Vf2Wizard : public TR {
     typedef TR Base;
     typedef typename TR::Graph1 Graph1;
@@ -507,12 +466,16 @@
   public:
     ///Constructor
-    Vf2Wizard(const Graph1 &g1,const Graph2 &g2) : Base(g1,g2) {}
+    Vf2Wizard(const Graph1 &g1,const Graph2 &g2) : Base(g1,g2) {
+    }
 
     ///Copy constructor
-    Vf2Wizard(const Base &b) : Base(b) {}
+    Vf2Wizard(const Base &b) : Base(b) { }
+
+    ///Copy constructor
+    Vf2Wizard(const Vf2Wizard &b) : Base(b) {}
 
 
     template<class T>
-    struct SetMappingBase : public Base {
+    struct SetMappingBase : public Base{
       typedef T Mapping;
       SetMappingBase(const Base &b) : Base(b) {}
@@ -525,6 +488,5 @@
     ///the map that stores the found embedding.
     template<class T>
-    Vf2Wizard< SetMappingBase<T> > mapping(const T &t)
-    {
+    Vf2Wizard< SetMappingBase<T> > mapping(const T &t) {
       Base::_mapping=reinterpret_cast<void*>(const_cast<T*>(&t));
       Base::_local_mapping = false;
@@ -537,5 +499,6 @@
       NodeEq _node_eq;
       SetNodeEqBase(const Base &b, const NE &node_eq)
-        : Base(b), _node_eq(node_eq) {}
+        : Base(b), _node_eq(node_eq){
+      }
     };
 
@@ -550,6 +513,5 @@
     ///always true operator.
     template<class T>
-    Vf2Wizard< SetNodeEqBase<T> > nodeEq(const T &node_eq)
-    {
+    Vf2Wizard< SetNodeEqBase<T> > nodeEq(const T &node_eq) {
       return Vf2Wizard<SetNodeEqBase<T> >(SetNodeEqBase<T>(*this,node_eq));
     }
@@ -561,7 +523,7 @@
     ///the node labels defining equivalence relation between them.
     ///
-    ///\param m1 It is arbitrary \ref concepts::ReadMap "readable node map"
+    ///\param m1 An arbitrary \ref concepts::ReadMap "readable node map"
     ///of g1.
-    ///\param m2 It is arbitrary \ref concepts::ReadMap "readable node map"
+    ///\param m2 An arbitrary \ref concepts::ReadMap "readable node map"
     ///of g2.
     ///
@@ -569,6 +531,5 @@
     template<class M1, class M2>
     Vf2Wizard< SetNodeEqBase<bits::vf2::MapEq<M1,M2> > >
-    nodeLabels(const M1 &m1,const M2 &m2)
-    {
+    nodeLabels(const M1 &m1,const M2 &m2){
       return nodeEq(bits::vf2::MapEq<M1,M2>(m1,m2));
     }
@@ -582,7 +543,6 @@
     ///The mapping type is set to \ref SUBGRAPH by default.
     ///
-    ///\sa See \ref Vf2MappingType for the possible values.
-    Vf2Wizard<Base> &mappingType(Vf2MappingType m_type)
-    {
+    ///\sa See \ref MappingType for the possible values.
+    Vf2Wizard<Base> &mappingType(MappingType m_type) {
       _mapping_type = m_type;
       return *this;
@@ -594,6 +554,5 @@
     ///\ref named-templ-param "Named parameter" for setting
     ///the mapping type to \ref INDUCED.
-    Vf2Wizard<Base> &induced()
-    {
+    Vf2Wizard<Base> &induced() {
       _mapping_type = INDUCED;
       return *this;
@@ -605,10 +564,10 @@
     ///\ref named-templ-param "Named parameter" for setting
     ///the mapping type to \ref ISOMORPH.
-    Vf2Wizard<Base> &iso()
-    {
+    Vf2Wizard<Base> &iso() {
       _mapping_type = ISOMORPH;
       return *this;
     }
 
+
     ///Runs VF2 algorithm.
 
@@ -616,7 +575,6 @@
     ///
     ///\retval true if a mapping is found.
-    ///\retval false if there is no (more) mapping.
-    bool run()
-    {
+    ///\retval false if there is no mapping.
+    bool run(){
       if(Base::_local_mapping)
         Base::createMapping();
@@ -631,4 +589,44 @@
       if(Base::_local_mapping)
         delete reinterpret_cast<Mapping*>(_mapping);
+
+      return ret;
+    }
+
+    ///Get a pointer to the generated Vf2 object.
+
+    ///Gives a pointer to the generated Vf2 object.
+    ///
+    ///\return Pointer to the generated Vf2 object.
+    ///\warning Don't forget to delete the referred Vf2 object after use.
+    Vf2<Graph1, Graph2, Mapping, NodeEq >* getPtrToVf2Object() {
+      if(Base::_local_mapping)
+        Base::createMapping();
+      Vf2<Graph1, Graph2, Mapping, NodeEq >* ptr =
+        new Vf2<Graph1, Graph2, Mapping, NodeEq>
+        (_g1, _g2, *reinterpret_cast<Mapping*>(_mapping), _node_eq);
+      ptr->mappingType(_mapping_type);
+      if(Base::_local_mapping)
+        ptr->_deallocMappingAfterUse = true;
+      return ptr;
+    }
+
+    ///Counts the number of mappings.
+
+    ///This method counts the number of mappings.
+    ///
+    /// \return The number of mappings.
+    int count() {
+      if(Base::_local_mapping)
+        Base::createMapping();
+
+      Vf2<Graph1, Graph2, Mapping, NodeEq>
+        alg(_g1, _g2, *reinterpret_cast<Mapping*>(_mapping), _node_eq);
+      if(Base::_local_mapping)
+        alg._deallocMappingAfterUse = true;
+      alg.mappingType(_mapping_type);
+
+      int ret = 0;
+      while(alg.find())
+        ++ret;
 
       return ret;
@@ -645,18 +643,30 @@
   ///The following examples show how to use these parameters.
   ///\code
-  ///  // Find an embedding of graph g into graph h
+  ///  // Find an embedding of graph g1 into graph g2
   ///  ListGraph::NodeMap<ListGraph::Node> m(g);
-  ///  vf2(g,h).mapping(m).run();
-  ///
-  ///  // Check whether graphs g and h are isomorphic
-  ///  bool is_iso = vf2(g,h).iso().run();
+  ///  vf2(g1,g2).mapping(m).run();
+  ///
+  ///  // Check whether graphs g1 and g2 are isomorphic
+  ///  bool is_iso = vf2(g1,g2).iso().run();
+  ///
+  ///  // Count the number of isomorphisms
+  ///  int num_isos = vf2(g1,g2).iso().count();
+  ///
+  ///  // Iterate through all the induced subgraph mappings of graph g1 into g2
+  ///  auto* myVf2 = vf2(g1,g2).mapping(m).nodeLabels(c1,c2)
+  ///  .induced().getPtrToVf2Object();
+  ///  while(myVf2->find()){
+  ///    //process the current mapping m
+  ///  }
+  ///  delete myVf22;
   ///\endcode
-  ///\warning Don't forget to put the \ref Vf2Wizard::run() "run()"
+  ///\warning Don't forget to put the \ref Vf2Wizard::run() "run()",
+  ///\ref Vf2Wizard::count() "count()" or
+  ///the \ref Vf2Wizard::getPtrToVf2Object() "getPtrToVf2Object()"
   ///to the end of the expression.
   ///\sa Vf2Wizard
   ///\sa Vf2
   template<class G1, class G2>
-  Vf2Wizard<Vf2WizardBase<G1,G2> > vf2(const G1 &g1, const G2 &g2)
-  {
+  Vf2Wizard<Vf2WizardBase<G1,G2> > vf2(const G1 &g1, const G2 &g2) {
     return Vf2Wizard<Vf2WizardBase<G1,G2> >(g1,g2);
   }
