Main MRPT website > C++ reference for MRPT 1.9.9
CGraphPartitioner_impl.h
Go to the documentation of this file.
1 /* +------------------------------------------------------------------------+
2  | Mobile Robot Programming Toolkit (MRPT) |
3  | http://www.mrpt.org/ |
4  | |
5  | Copyright (c) 2005-2018, Individual contributors, see AUTHORS file |
6  | See: http://www.mrpt.org/Authors - All rights reserved. |
7  | Released under BSD License. See details in http://www.mrpt.org/License |
8  +------------------------------------------------------------------------+ */
9 
10 #if !defined(CGRAPHPARTITIONER_H)
11 #error "This file can't be included from outside of CGraphPartitioner.h"
12 #endif
13 
14 namespace mrpt
15 {
16 namespace graphs
17 {
18 /*---------------------------------------------------------------
19  SpectralPartition
20  ---------------------------------------------------------------*/
21 template <class GRAPH_MATRIX, typename num_t>
22 void CGraphPartitioner<GRAPH_MATRIX, num_t>::SpectralBisection(
23  GRAPH_MATRIX& in_A, std::vector<uint32_t>& out_part1,
24  std::vector<uint32_t>& out_part2, num_t& out_cut_value, bool forceSimetry)
25 {
26  size_t nodeCount; // Nodes count
27  GRAPH_MATRIX Adj, eigenVectors, eigenValues;
28 
29  // Check matrix is square:
30  if (in_A.cols() != int(nodeCount = in_A.rows()))
31  THROW_EXCEPTION("Weights matrix is not square!!");
32 
33  // Shi & Malik's method
34  //-----------------------------------------------------------------
35 
36  // forceSimetry?
37  if (forceSimetry)
38  {
39  Adj.setSize(nodeCount, nodeCount);
40  for (size_t i = 0; i < nodeCount; i++)
41  for (size_t j = i; j < nodeCount; j++)
42  Adj(i, j) = Adj(j, i) = 0.5f * (in_A(i, j) + in_A(j, i));
43  }
44  else
45  Adj = in_A;
46 
47  // Compute eigen-vectors of laplacian:
48  GRAPH_MATRIX LAPLACIAN;
49  Adj.laplacian(LAPLACIAN);
50 
51  LAPLACIAN.eigenVectors(eigenVectors, eigenValues);
52 
53  // Execute the bisection
54  // ------------------------------------
55  // Second smallest eigen-vector
56  double mean = 0;
57  size_t colNo = 1; // second smallest
58  size_t nRows = eigenVectors.rows();
59 
60  // for (i=0;i<eigenVectors.cols();i++) mean+=eigenVectors(colNo,i);
61  for (size_t i = 0; i < nRows; i++) mean += eigenVectors(i, colNo);
62  mean /= nRows;
63 
64  out_part1.clear();
65  out_part2.clear();
66 
67  for (size_t i = 0; i < nRows; i++)
68  {
69  if (eigenVectors(i, colNo) >= mean)
70  out_part1.push_back(i);
71  else
72  out_part2.push_back(i);
73  }
74 
75  // Special and strange case: Constant eigenvector: Split nodes in two
76  // equally sized parts arbitrarily:
77  // ------------------------------------------------------------------
78  if (!out_part1.size() || !out_part2.size())
79  {
80  out_part1.clear();
81  out_part2.clear();
82  // Assign 50%-50%:
83  for (int i = 0; i < Adj.cols(); i++)
84  if (i <= Adj.cols() / 2)
85  out_part1.push_back(i);
86  else
87  out_part2.push_back(i);
88  }
89 
90  // Compute the N-cut value
91  out_cut_value = nCut(Adj, out_part1, out_part2);
92 }
93 
94 /*---------------------------------------------------------------
95  RecursiveSpectralPartition
96  ---------------------------------------------------------------*/
97 template <class GRAPH_MATRIX, typename num_t>
98 void CGraphPartitioner<GRAPH_MATRIX, num_t>::RecursiveSpectralPartition(
99  GRAPH_MATRIX& in_A, std::vector<std::vector<uint32_t>>& out_parts,
100  num_t threshold_Ncut, bool forceSimetry, bool useSpectralBisection,
101  bool recursive, unsigned minSizeClusters, const bool verbose)
102 {
103  MRPT_START
104 
105  size_t nodeCount;
106  std::vector<uint32_t> p1, p2;
107  num_t cut_value;
108  size_t i, j;
109  GRAPH_MATRIX Adj;
110 
111  out_parts.clear();
112 
113  // Check matrix is square:
114  if (in_A.cols() != int(nodeCount = in_A.rows()))
115  THROW_EXCEPTION("Weights matrix is not square!!");
116 
117  if (nodeCount == 1)
118  {
119  // Don't split, there is just a node!
120  p1.push_back(0);
121  out_parts.push_back(p1);
122  return;
123  }
124 
125  // forceSimetry?
126  if (forceSimetry)
127  {
128  Adj.setSize(nodeCount, nodeCount);
129  for (i = 0; i < nodeCount; i++)
130  for (j = i; j < nodeCount; j++)
131  Adj(i, j) = Adj(j, i) = 0.5f * (in_A(i, j) + in_A(j, i));
132  }
133  else
134  Adj = in_A;
135 
136  // Make bisection
137  if (useSpectralBisection)
138  SpectralBisection(Adj, p1, p2, cut_value, false);
139  else
140  exactBisection(Adj, p1, p2, cut_value, false);
141 
142  if (verbose)
143  std::cout << format(
144  "Cut:%u=%u+%u,nCut=%.02f->", (unsigned int)nodeCount,
145  (unsigned int)p1.size(), (unsigned int)p2.size(), cut_value);
146 
147  // Is it a useful partition?
148  if (cut_value > threshold_Ncut || p1.size() < minSizeClusters ||
149  p2.size() < minSizeClusters)
150  {
151  if (verbose) std::cout << "->NO!" << std::endl;
152 
153  // No:
154  p1.clear();
155  for (i = 0; i < nodeCount; i++) p1.push_back(i);
156  out_parts.push_back(p1);
157  }
158  else
159  {
160  if (verbose) std::cout << "->YES!" << std::endl;
161 
162  // Yes:
163  std::vector<std::vector<uint32_t>> p1_parts, p2_parts;
164 
165  if (recursive)
166  {
167  // Split "p1":
168  // --------------------------------------------
169  // sub-matrix:
170  GRAPH_MATRIX A_1(p1.size(), p1.size());
171  for (i = 0; i < p1.size(); i++)
172  for (j = 0; j < p1.size(); j++) A_1(i, j) = in_A(p1[i], p1[j]);
173 
174  RecursiveSpectralPartition(
175  A_1, p1_parts, threshold_Ncut, forceSimetry,
176  useSpectralBisection, recursive, minSizeClusters);
177 
178  // Split "p2":
179  // --------------------------------------------
180  // sub-matrix:
181  GRAPH_MATRIX A_2(p2.size(), p2.size());
182  for (i = 0; i < p2.size(); i++)
183  for (j = 0; j < p2.size(); j++) A_2(i, j) = in_A(p2[i], p2[j]);
184 
185  RecursiveSpectralPartition(
186  A_2, p2_parts, threshold_Ncut, forceSimetry,
187  useSpectralBisection, recursive, minSizeClusters);
188 
189  // Build "out_parts" from "p1_parts" + "p2_parts"
190  // taken care of indexes mapping!
191  // -------------------------------------------------------------------------------------
192  // remap p1 nodes:
193  for (i = 0; i < p1_parts.size(); i++)
194  {
195  for (j = 0; j < p1_parts[i].size(); j++)
196  p1_parts[i][j] = p1[p1_parts[i][j]];
197  out_parts.push_back(p1_parts[i]);
198  }
199 
200  // remap p2 nodes:
201  for (i = 0; i < p2_parts.size(); i++)
202  {
203  for (j = 0; j < p2_parts[i].size(); j++)
204  p2_parts[i][j] = p2[p2_parts[i][j]];
205 
206  out_parts.push_back(p2_parts[i]);
207  }
208  }
209  else
210  {
211  // Force bisection only:
212  out_parts.clear();
213  out_parts.push_back(p1);
214  out_parts.push_back(p2);
215  }
216  }
217 
218  MRPT_END
219 }
220 
221 /*---------------------------------------------------------------
222  nCut
223  ---------------------------------------------------------------*/
224 template <class GRAPH_MATRIX, typename num_t>
225 num_t CGraphPartitioner<GRAPH_MATRIX, num_t>::nCut(
226  const GRAPH_MATRIX& in_A, const std::vector<uint32_t>& in_part1,
227  const std::vector<uint32_t>& in_part2)
228 {
229  unsigned int i, j;
230  size_t size1 = in_part1.size();
231  size_t size2 = in_part2.size();
232 
233  // Compute the N-cut value
234  // -----------------------------------------------
235  num_t cut_AB = 0;
236  for (i = 0; i < size1; i++)
237  for (j = 0; j < size2; j++) cut_AB += in_A(in_part1[i], in_part2[j]);
238 
239  num_t assoc_AA = 0;
240  for (i = 0; i < size1; i++)
241  for (j = i; j < size1; j++)
242  if (i != j) assoc_AA += in_A(in_part1[i], in_part1[j]);
243 
244  num_t assoc_BB = 0;
245  for (i = 0; i < size2; i++)
246  for (j = i; j < size2; j++)
247  if (i != j) assoc_BB += in_A(in_part2[i], in_part2[j]);
248 
249  num_t assoc_AV = assoc_AA + cut_AB;
250  num_t assoc_BV = assoc_BB + cut_AB;
251 
252  if (!cut_AB)
253  return 0;
254  else
255  return cut_AB / assoc_AV + cut_AB / assoc_BV;
256 }
257 
258 /*---------------------------------------------------------------
259  exactBisection
260  ---------------------------------------------------------------*/
261 template <class GRAPH_MATRIX, typename num_t>
262 void CGraphPartitioner<GRAPH_MATRIX, num_t>::exactBisection(
263  GRAPH_MATRIX& in_A, std::vector<uint32_t>& out_part1,
264  std::vector<uint32_t>& out_part2, num_t& out_cut_value, bool forceSimetry)
265 {
266  size_t nodeCount; // Nodes count
267  size_t i, j;
268  GRAPH_MATRIX Adj;
269  std::vector<bool> partition, bestPartition;
270  std::vector<uint32_t> part1, part2;
271  num_t partCutValue, bestCutValue = std::numeric_limits<num_t>::max();
272  bool end = false;
273 
274  // Check matrix is square:
275  if (in_A.cols() != int(nodeCount = in_A.rows()))
276  THROW_EXCEPTION("Weights matrix is not square!!");
277 
278  ASSERT_(nodeCount >= 2);
279 
280  // forceSimetry?
281  if (forceSimetry)
282  {
283  Adj.setSize(nodeCount, nodeCount);
284  for (i = 0; i < nodeCount; i++)
285  for (j = i; j < nodeCount; j++)
286  Adj(i, j) = Adj(j, i) = 0.5f * (in_A(i, j) + in_A(j, i));
287  }
288  else
289  Adj = in_A;
290 
291  // Brute force: compute all posible partitions:
292  //-----------------------------------------------------------------
293 
294  // First combination: 1000...0000
295  // Last combination: 1111...1110
296  partition.clear();
297  partition.resize(nodeCount, false);
298  partition[0] = true;
299 
300  while (!end)
301  {
302  // Build partitions from binary vector:
303  part1.clear();
304  part2.clear();
305 
306  for (i = 0; i < nodeCount; i++)
307  {
308  if (partition[i])
309  part2.push_back(i);
310  else
311  part1.push_back(i);
312  }
313 
314  // Compute the n-cut:
315  partCutValue = nCut(Adj, part1, part2);
316 
317  if (partCutValue < bestCutValue)
318  {
319  bestCutValue = partCutValue;
320  bestPartition = partition;
321  }
322 
323  // Next combo in the binary vector:
324  i = 0;
325  bool carry = false;
326  do
327  {
328  carry = partition[i];
329  partition[i] = !partition[i];
330  i++;
331  } while (carry && i < nodeCount);
332 
333  // End criterion:
334  end = true;
335  for (i = 0; end && i < nodeCount; i++) end = end && partition[i];
336 
337  }; // End of while
338 
339  // Return the best partition:
340  out_cut_value = bestCutValue;
341 
342  out_part1.clear();
343  out_part2.clear();
344 
345  for (i = 0; i < nodeCount; i++)
346  {
347  if (bestPartition[i])
348  out_part2.push_back(i);
349  else
350  out_part1.push_back(i);
351  }
352 }
353 
354 } // namespace graphs
355 } // namespace mrpt
end
GLuint GLuint end
Definition: glext.h:3528
mrpt
This is the global namespace for all Mobile Robot Programming Toolkit (MRPT) libraries.
Definition: CKalmanFilterCapable.h:30
THROW_EXCEPTION
#define THROW_EXCEPTION(msg)
Definition: exceptions.h:41
ASSERT_
#define ASSERT_(f)
Defines an assertion mechanism.
Definition: exceptions.h:113
eigenVectors
EIGEN_STRONG_INLINE bool eigenVectors(MATRIX1 &eVecs, MATRIX2 &eVals) const
[For square matrices only] Compute the eigenvectors and eigenvalues (sorted), both returned as matric...
mrpt::graphs::CGraphPartitioner::nCut
static num_t nCut(const GRAPH_MATRIX &in_A, const std::vector< uint32_t > &in_part1, const std::vector< uint32_t > &in_part2)
Returns the normaliced cut of a graph, given its adjacency matrix A and a bisection:
Definition: CGraphPartitioner_impl.h:225
MRPT_START
#define MRPT_START
Definition: exceptions.h:262
mrpt::format
std::string format(const char *fmt,...) MRPT_printf_format_check(1
A std::string version of C sprintf.
Definition: format.cpp:16
eigenValues
EIGEN_STRONG_INLINE void eigenValues(VECTOR &eVals) const
[For square matrices only] Compute the eigenvectors and eigenvalues (sorted), and return only the eig...
Definition: eigen_plugins.h:868
mrpt::graphs::CGraphPartitioner::RecursiveSpectralPartition
static void RecursiveSpectralPartition(GRAPH_MATRIX &in_A, std::vector< std::vector< uint32_t >> &out_parts, num_t threshold_Ncut=1, bool forceSimetry=true, bool useSpectralBisection=true, bool recursive=true, unsigned minSizeClusters=1, const bool verbose=false)
Performs the spectral recursive partition into K-parts for a given graph.
Definition: CGraphPartitioner_impl.h:98
mean
EIGEN_STRONG_INLINE double mean() const
Computes the mean of the entire matrix.
Definition: eigen_plugins.h:427
MRPT_END
#define MRPT_END
Definition: exceptions.h:266
mrpt::graphs::CGraphPartitioner::SpectralBisection
static void SpectralBisection(GRAPH_MATRIX &in_A, std::vector< uint32_t > &out_part1, std::vector< uint32_t > &out_part2, num_t &out_cut_value, bool forceSimetry=true)
Performs the spectral bisection of a graph.
Definition: CGraphPartitioner_impl.h:22
mrpt::graphs::CGraphPartitioner::exactBisection
static void exactBisection(GRAPH_MATRIX &in_A, std::vector< uint32_t > &out_part1, std::vector< uint32_t > &out_part2, num_t &out_cut_value, bool forceSimetry=true)
Performs an EXACT minimum n-Cut graph bisection, (Use CGraphPartitioner::SpectralBisection for a fast...
Definition: CGraphPartitioner_impl.h:262


Page generated by Doxygen 1.8.17 for MRPT 1.9.9 Git: ad3a9d8ae Tue May 1 23:10:22 2018 -0700 at mié 12 jul 2023 10:03:34 CEST