inference.blockmodel: simplify and speedup vertex movements

src/graph/inference/graph_blockmodel.cc

@@ -153,6 +153,8 @@ void export_blockmodel_state()
                  &state_t::remove_vertex;
              void (state_t::*add_vertex)(size_t, size_t) =
                  &state_t::add_vertex;
+             void (state_t::*move_vertex)(size_t, size_t) =
+                 &state_t::move_vertex;
              void (state_t::*remove_vertices)(python::object) =
                  &state_t::remove_vertices;
              void (state_t::*add_vertices)(python::object, python::object) =
@@ -180,7 +182,7 @@ void export_blockmodel_state()
                  .def("add_vertex", add_vertex)
                  .def("remove_vertices", remove_vertices)
                  .def("add_vertices", add_vertices)
-                 .def("move_vertex", &state_t::move_vertex)
+                 .def("move_vertex", move_vertex)
                  .def("move_vertices", move_vertices)
                  .def("set_partition", set_partition)
                  .def("virtual_move", virtual_move)

src/graph/inference/graph_blockmodel_overlap.hh

@@ -86,7 +86,7 @@ public:
           _c_mrs(_mrs.get_checked()),
           _c_brec(_brec.get_checked()),
           _c_bdrec(_bdrec.get_checked()),
-          _emat(_g, _b, _bg, rng),
+          _emat(_bg, rng),
           _overlap_stats(_g, _b, _half_edges, _node_index, num_vertices(_bg))
     {
         for (auto r : vertices_range(_bg))
@@ -117,12 +117,13 @@ public:
         {
             size_t s = _b[u];
 
-            auto e = *out_edges(v, _g).first;
-            auto& me = _emat.get_bedge(e);
+            auto& me = _emat.get_me(r, s);
 
             _mrs[me] -= 1;
             _mrp[r] -= 1;
             _mrm[s] -= 1;
+
+            auto e = *out_edges(v, _g).first;
             eop(e, me);
 
             assert(_mrs[me] >= 0);
@@ -135,12 +136,13 @@ public:
         {
             size_t s = _b[u];
 
-            auto e = *in_edge_iteratorS<g_t>().get_edges(v, _g).first;
-            auto& me = _emat.get_bedge(e);
+            auto& me = _emat.get_me(s, r);
 
             _mrs[me] -= 1;
             _mrp[s] -= 1;
             _mrm[r] -= 1;
+
+            auto e = *in_edge_iteratorS<g_t>().get_edges(v, _g).first;
             eop(e, me);
 
             if (_mrs[me] == 0)
@@ -199,14 +201,13 @@ public:
                 _c_bdrec[me] = 0;
             }
 
-            auto e = *out_edges(v, _g).first;
-            _emat.get_bedge(e) = me;
-
             assert(me == _emat.get_me(r, s));
 
             _mrs[me] += 1;
             _mrp[r] += 1;
             _mrm[s] += 1;
+
+            auto e = *out_edges(v, _g).first;
             eop(e, me);
         }
 
@@ -226,14 +227,12 @@ public:
                 _c_bdrec[me] = 0;
             }
 
-            auto e = *in_edge_iteratorS<g_t>().get_edges(v, _g).first;
-            _emat.get_bedge(e) = me;
-
             assert(me == _emat.get_me(s, r));
 
             _mrs[me] += 1;
             _mrp[s] += 1;
             _mrm[r] += 1;
+            auto e = *in_edge_iteratorS<g_t>().get_edges(v, _g).first;
             eop(e, me);
         }
 
@@ -331,9 +330,10 @@ public:
     {
         auto mv_entries = [&](auto&&... args)
             {
-                move_entries(v, r, nr, _b, _emat.get_bedge_map(), _g, _bg,
-                             m_entries, is_loop_overlap(_overlap_stats),
-                             _eweight, args...);
+                move_entries(v, r, nr, [&](auto u) { return this->_b[u]; },
+                             _g, _eweight, m_entries,
+                             [](auto) { return false; },
+                             is_loop_overlap(_overlap_stats), args...);
             };
 
         switch (_rec_type)
@@ -359,9 +359,6 @@ public:
         if (r == nr)
             return 0.;
 
-        m_entries.clear();
-        get_move_entries(v, r, nr, m_entries);
-
         size_t kout = out_degreeS()(v, _g);
         size_t kin = 0;
         if (is_directed::apply<g_t>::type::value)
@@ -424,6 +421,9 @@ public:
         if (!allow_move(r, nr))
             return std::numeric_limits<double>::infinity();
 
+        m_entries.clear();
+        get_move_entries(v, r, nr, m_entries);
+
         double dS = 0;
         if (ea.adjacency)
         {
@@ -432,11 +432,6 @@ public:
             else
                 dS = virtual_move_sparse<false>(v, nr, ea.multigraph, m_entries);
         }
-        else
-        {
-            m_entries.clear();
-            get_move_entries(v, r, nr, m_entries);
-        }
 
         if (ea.partition_dl || ea.degree_dl || ea.edges_dl)
         {
@@ -453,7 +448,7 @@ public:
         switch (_rec_type)
         {
         case weight_type::POSITIVE: // positive weights
-            entries_op(m_entries,
+            entries_op(m_entries, _emat,
                        [&](auto, auto, auto& me, auto& delta)
                        {
                            size_t ers = 0;
@@ -472,7 +467,7 @@ public:
                        });
             break;
         case weight_type::SIGNED: // positive and negative weights
-            entries_op(m_entries,
+            entries_op(m_entries, _emat,
                        [&](auto, auto, auto& me, auto& delta)
                        {
                            size_t ers = 0;
@@ -613,18 +608,20 @@ public:
                 size_t ew = _eweight[e];
                 w += ew;
 
-                int mts;
-                if (t == r && s == size_t(_b[u]))
-                    mts = _mrs[_emat.get_bedge(e)];
-                else
-                    mts = get_beprop(t, s, _mrs, _emat);
+                int mts = 0;
+                const auto& me = m_entries.get_me(t, s, _emat);
+                if (me != _emat.get_null_edge())
+                    mts = _mrs[me];
                 int mtp = _mrp[t];
                 int mst = mts;
                 int mtm = mtp;
 
                 if (is_directed::apply<g_t>::type::value)
                 {
-                    mst = get_beprop(s, t, _mrs, _emat);
+                    mst = 0;
+                    const auto& me = m_entries.get_me(s, t, _emat);
+                    if (me != _emat.get_null_edge())
+                        mst = _mrs[me];
                     mtm = _mrm[t];
                 }
 
@@ -1005,8 +1002,8 @@ public:
     typename bdrec_t::checked_t _c_bdrec;
 
     typedef typename std::conditional<use_hash_t::value,
-                                      EHash<g_t, bg_t>,
-                                      EMat<g_t, bg_t>>::type
+                                      EHash<bg_t>,
+                                      EMat<bg_t>>::type
         emat_t;
     emat_t _emat;
 

src/graph/inference/graph_blockmodel_overlap_util.hh

@@ -1236,19 +1236,18 @@ class SingleEntrySet
 public:
     typedef typename graph_traits<BGraph>::edge_descriptor bedge_t;
 
-    SingleEntrySet() : _pos(0) {}
+    SingleEntrySet() : _pos(0), _mes_pos(0) {}
     SingleEntrySet(size_t) : SingleEntrySet() {}
 
     void set_move(size_t, size_t) {}
 
     template <class... DVals>
-    void insert_delta(size_t t, size_t s, const bedge_t& me, DVals... delta)
+    void insert_delta(size_t t, size_t s, DVals... delta)
     {
         if (!is_directed::apply<Graph>::type::value && (t > s))
             std::swap(t, s);
         _entries[_pos] = make_pair(t, s);
         add_to_tuple(_delta[_pos], delta...);
-        _mes[_pos] = me;
         ++_pos;
     }
 
@@ -1270,18 +1269,51 @@ public:
         for (auto& d : _delta)
             d = std::tuple<EVals...>();
         _pos = 0;
+        _mes_pos = 0;
     }
 
     const std::array<pair<size_t, size_t>,2>& get_entries() const { return _entries; }
     const std::array<std::tuple<EVals...>, 2>& get_delta() const { return _delta; }
-    std::array<bedge_t, 2>& get_mes() { return _mes; }
     const bedge_t& get_null_edge() const { return _null_edge; }
 
+    template <class Emat>
+    std::array<bedge_t, 2>& get_mes(Emat& emat)
+    {
+        for (; _mes_pos < 2; ++_mes_pos)
+        {
+            auto& entry = _entries[_mes_pos];
+            _mes[_mes_pos] = emat.get_me(entry.first, entry.second);
+        }
+        return _mes;
+    }
+
+    template <class Emat>
+    const bedge_t& get_me(size_t t, size_t s, Emat& emat)
+    {
+        if (!is_directed::apply<Graph>::type::value && (t > s))
+            std::swap(t, s);
+        for (size_t i = 0; i < 2; ++i)
+        {
+            auto& entry = _entries[i];
+            if (entry.first == t && entry.second == s)
+            {
+                if (i >= _mes_pos)
+                {
+                    _mes[i] = emat.get_me(t, s);
+                    _mes_pos++;
+                }
+                return _mes[i];
+            }
+        }
+        return emat.get_me(t, s);
+    }
+
 private:
     size_t _pos;
     std::array<pair<size_t, size_t>, 2> _entries;
     std::array<std::tuple<EVals...>, 2> _delta;
     std::array<bedge_t, 2> _mes;
+    size_t _mes_pos;
 
     static const std::tuple<EVals...> _null_delta;
     static const bedge_t _null_edge;

src/graph/inference/graph_blockmodel_util.hh

@@ -768,19 +768,17 @@ private:
 // we're using an adjacency list to store the block structure (it is simply an
 // adjacency matrix)
 
-template <class Graph, class BGraph>
+template <class BGraph>
 class EMat
 {
 public:
-    template <class Vprop, class RNG>
-    EMat(Graph& g, Vprop b, BGraph& bg, RNG&)
-        : _bedge(get(edge_index_t(), g), 0)
+    template <class RNG>
+    EMat(BGraph& bg, RNG&)
     {
-        sync(g, b, bg);
+        sync(bg);
     }
 
-    template <class Vprop>
-    void sync(Graph& g, Vprop b, BGraph& bg)
+    void sync(BGraph& bg)
     {
         size_t B = num_vertices(bg);
         _mat.resize(boost::extents[B][B]);
@@ -793,15 +791,6 @@ public:
             if (!is_directed::apply<BGraph>::type::value)
                 _mat[target(e, bg)][source(e, bg)] = e;
         }
-
-        auto bedge_c = _bedge.get_checked();
-        for (auto e : edges_range(g))
-        {
-            auto r = b[source(e, g)];
-            auto s = b[target(e, g)];
-            bedge_c[e] = _mat[r][s];
-            assert(bedge_c[e] != _null_edge);
-        }
     }
 
     typedef typename graph_traits<BGraph>::vertex_descriptor vertex_t;
@@ -825,30 +814,21 @@ public:
         if (delete_edge)
         {
             _mat[r][s] = _null_edge;
-            if (!is_directed::apply<Graph>::type::value)
+            if (!is_directed::apply<BGraph>::type::value)
                 _mat[s][r] = _null_edge;
             remove_edge(me, bg);
         }
     }
 
-    template <class Edge>
-    auto& get_bedge(const Edge& e) { return _bedge[e]; }
-    auto& get_bedge_map() { return _bedge; }
-    const auto& get_bedge_map() const { return _bedge; }
     const auto& get_null_edge() const { return _null_edge; }
 
 private:
     multi_array<edge_t, 2> _mat;
-    typedef typename property_map_type::apply
-        <edge_t,
-         typename property_map<Graph,
-                               edge_index_t>::type>::type bedge_t;
-    typename bedge_t::unchecked_t _bedge;
     static const edge_t _null_edge;
 };
 
-template <class Graph, class BGraph>
-const typename EMat<Graph, BGraph>::edge_t EMat<Graph, BGraph>::_null_edge;
+template <class BGraph>
+const typename EMat<BGraph>::edge_t EMat<BGraph>::_null_edge;
 
 
 template <class Key>
@@ -875,22 +855,20 @@ private:
 // we're using an adjacency list to store the block structure (this is like
 // EMat above, but takes less space and is slower)
 
-template <class Graph, class BGraph>
+template <class BGraph>
 class EHash
 {
 public:
 
-    template <class Vprop, class RNG>
-    EHash(Graph& g, Vprop b, BGraph& bg, RNG& rng)
+    template <class RNG>
+    EHash(BGraph& bg, RNG& rng)
         : _hash_function(num_vertices(bg), rng),
-          _hash(num_vertices(bg), ehash_t(0, _hash_function)),
-          _bedge(get(edge_index_t(), g), 0)
+          _hash(num_vertices(bg), ehash_t(0, _hash_function))
     {
-        sync(g, b, bg);
+        sync(bg);
     }
 
-    template <class Vprop>
-    void sync(Graph& g, Vprop b, BGraph& bg)
+    void sync(BGraph& bg)
     {
         _hash.clear();
         _hash.resize(num_vertices(bg), ehash_t(0, _hash_function));
@@ -900,15 +878,6 @@ public:
             assert(get_me(source(e, bg), target(e, bg)) == _null_edge);
             put_me(source(e, bg), target(e, bg), e);
         }
-
-        auto bedge_c = _bedge.get_checked();
-        for (auto e : edges_range(g))
-        {
-            auto r = b[source(e, g)];
-            auto s = b[target(e, g)];
-            bedge_c[e] = get_me(r, s);
-            assert(bedge_c[e] != _null_edge);
-        }
     }
 
     typedef typename graph_traits<BGraph>::vertex_descriptor vertex_t;
@@ -916,6 +885,8 @@ public:
 
     const auto& get_me(vertex_t r, vertex_t s) const
     {
+        if (!is_directed::apply<BGraph>::type::value && r > s)
+            std::swap(r, s);
         auto& map = _hash[r];
         const auto& iter = map.find(s);
         if (iter == map.end())
@@ -925,10 +896,10 @@ public:
 
     void put_me(vertex_t r, vertex_t s, const edge_t& e)
     {
+        if (!is_directed::apply<BGraph>::type::value && r > s)
+            std::swap(r, s);
         assert(r < _hash.size());
         _hash[r][s] = e;
-        if (!is_directed::apply<Graph>::type::value)
-            _hash[s][r] = e;
     }
 
     void remove_me(vertex_t r, vertex_t s, const edge_t& me, BGraph& bg,
@@ -936,34 +907,25 @@ public:
     {
         if (delete_edge)
         {
+            if (!is_directed::apply<BGraph>::type::value && r > s)
+                std::swap(r, s);
             assert(r < _hash.size());
             _hash[r].erase(s);
-            if (!is_directed::apply<Graph>::type::value)