🥭 w/ louvainSeq, toleranceDeclineFactor

src/louvain.hxx

@@ -14,12 +14,13 @@ struct LouvainOptions {
   int repeat;
   T   resolution;
   T   tolerance;
-  T   phaseTolerance;
+  T   passTolerance;
+  T   toleranceDeclineFactor;
   int maxIterations;
-  int maxPhaseIterations;
+  int maxPasses;
 
-  LouvainOptions(int repeat=1, T resolution=1, T tolerance=0, T phaseTolerance=0, int maxIterations=500, int maxPhaseIterations=500) :
-  repeat(repeat), resolution(resolution), tolerance(tolerance), phaseTolerance(phaseTolerance), maxIterations(maxIterations), maxPhaseIterations(maxPhaseIterations) {}
+  LouvainOptions(int repeat=1, T resolution=1, T tolerance=1e-2, T passTolerance=0, T toleranceDeclineFactor=10, int maxIterations=500, int maxPasses=500) :
+  repeat(repeat), resolution(resolution), tolerance(tolerance), passTolerance(passTolerance), toleranceDeclineFactor(toleranceDeclineFactor), maxIterations(maxIterations), maxPasses(maxPasses) {}
 };
 
 

src/louvainOmp.hxx

@@ -160,7 +160,44 @@ void louvainChangeCommunity(vector<K>& vcom, vector<V>& ctot, const G& x, K u, K
  * @returns iterations
  */
 template <bool O, class G, class K, class V>
-int louvainMove(vector<K>& vdom, vector<K>& vcom, vector<V>& dtot, vector<V>& ctot, vector2d<K>& vcs, vector2d<V>& vcout, const G& x, const vector<V>& vtot, V M, V R, V E, int L) {
+int louvainMove(vector<K>& vdom, vector<K>& vcom, vector<V>& dtot, vector<V>& ctot, vector<K>& vcs, vector<V>& vcout, const G& x, const vector<V>& vtot, V M, V R, V E, int L) {
+  K S = x.span(), l = 0; V Q = V();
+  for (; l<L;) {
+    V el = V();
+    if (!O) copyValues(vcom, vdom);
+    if (!O) copyValues(ctot, dtot);
+    x.forEachVertexKey([&](auto u) {
+      louvainClearScan(vcs, vcout);
+      louvainScanCommunities(vcs, vcout, x, u, vdom);
+      auto [c, e] = louvainChooseCommunity(x, u, vdom, vtot, dtot, vcs, vcout, M, R);
+      if (O) { if (c)              louvainChangeCommunity(vcom, ctot, x, u, c, vdom, vtot); }
+      else   { if (c && c<vdom[u]) louvainChangeCommunity(vcom, ctot, x, u, c, vdom, vtot); }
+      el += e;   // l1-norm
+    }); ++l;
+    if (el<=E) break;
+  }
+  return l;
+}
+
+
+/**
+ * Louvain algorithm's local moving phase.
+ * @param vdom community each vertex belonged to (zeros)
+ * @param vcom community each vertex belongs to (initial)
+ * @param dtot previous total edge weight of each community (zeros)
+ * @param ctot total edge weight of each community (precalculated)
+ * @param vcs communities vertex u is linked to (temporary buffer)
+ * @param vcout total edge weight from vertex u to community C (temporary buffer)
+ * @param x original graph
+ * @param vtot total edge weight of each vertex
+ * @param M total weight of "undirected" graph (1/2 of directed graph)
+ * @param R resolution (0, 1]
+ * @param E tolerance (0)
+ * @param L max iterations (500)
+ * @returns iterations
+ */
+template <bool O, class G, class K, class V>
+int louvainMoveOmp(vector<K>& vdom, vector<K>& vcom, vector<V>& dtot, vector<V>& ctot, vector2d<K>& vcs, vector2d<V>& vcout, const G& x, const vector<V>& vtot, V M, V R, V E, int L) {
   K S = x.span(), l = 0; V Q = V();
   for (; l<L;) {
     V el = V();
@@ -226,13 +263,56 @@ void louvainLookupCommunities(vector<K>& a, const vector<K>& vcom) {
 }
 
 
+template <bool O, class G, class V=float>
+auto louvainSeq(const G& x, LouvainOptions<V> o={}) {
+  using K = typename G::key_type;
+  V   R = o.resolution;
+  V   D = o.passTolerance;
+  int L = o.maxIterations, l = 0;
+  int P = o.maxPasses, p = 0;
+  V   M = edgeWeight(x)/2;
+  size_t S = x.span();
+  vector<K> vdom(S), vcom(S), vcs, a(S);
+  vector<V> vtot(S), dtot(S), ctot(S), vcout(S);
+  float t = measureDurationMarked([&](auto mark) {
+    V E  = o.tolerance;
+    V Q0 = modularity(x, M, R);
+    G y  = duplicate(x);
+    fillValueU(vcom, K());
+    fillValueU(vtot, V());
+    fillValueU(ctot, V());
+    mark([&]() {
+      louvainVertexWeights(vtot, y);
+      louvainInitialize(vcom, ctot, y, vtot);
+      copyValues(vcom, a);
+      for (l=0, p=0; p<P;) {
+        if (O) l += louvainMove<O>(vcom, vcom, ctot, ctot, vcs, vcout, y, vtot, M, R, E, L);
+        else   l += louvainMove<O>(vdom, vcom, dtot, ctot, vcs, vcout, y, vtot, M, R, E, L);
+        y = louvainAggregate(y, vcom); ++p;
+        louvainLookupCommunities(a, vcom);
+        V Q = modularity(y, M, R);
+        V M = edgeWeight(y)/2;
+        if (Q-Q0<=D) break;
+        fillValueU(vcom, K());
+        fillValueU(vtot, V());
+        fillValueU(ctot, V());
+        louvainVertexWeights(vtot, y);
+        louvainInitialize(vcom, ctot, y, vtot);
+        E /= o.toleranceDeclineFactor;
+        Q0 = Q;
+      }
+    });
+  }, o.repeat);
+  return LouvainResult<K>(a, l, p, t);
+}
+
+
 template <bool O, class G, class V=float>
 auto louvainOmp(const G& x, LouvainOptions<V> o={}) {
   using K = typename G::key_type;
   V   R = o.resolution;
-  V   E = o.tolerance;
   V   D = o.passTolerance;
-  int L = o.maxIterations,      l = 0;
+  int L = o.maxIterations, l = 0;
   int P = o.maxPasses, p = 0;
   V   M = edgeWeight(x)/2;
   size_t S = x.span();
@@ -242,6 +322,7 @@ auto louvainOmp(const G& x, LouvainOptions<V> o={}) {
   vector2d<K> vcs(T);
   vector2d<V> vcout(T, vector<V>(S));
   float t = measureDurationMarked([&](auto mark) {
+    V E  = o.tolerance;
     V Q0 = modularity(x, M, R);
     G y  = duplicate(x);
     fillValueU(vcom, K());
@@ -264,6 +345,7 @@ auto louvainOmp(const G& x, LouvainOptions<V> o={}) {
         fillValueU(ctot, V());
         louvainVertexWeights(vtot, y);
         louvainInitialize(vcom, ctot, y, vtot);
+        E /= o.toleranceDeclineFactor;
         Q0 = Q;
       }
     });

src/main.hxx

@@ -14,4 +14,5 @@
 #include "properties.hxx"
 #include "modularity.hxx"
 #include "louvain.hxx"
+#include "louvainSeq.hxx"
 #include "louvainOmp.hxx"