COMBINATORIAL_BLAS 1.6
 
Loading...
Searching...
No Matches
BPMaximumMatching.cpp
Go to the documentation of this file.
1
2#ifdef THREADED
3#ifndef _OPENMP
4#define _OPENMP
5#endif
6
7#include <omp.h>
8int cblas_splits = 1;
9#endif
10
11#include "CombBLAS/CombBLAS.h"
12#include <mpi.h>
13#include <sys/time.h>
14#include <iostream>
15#include <functional>
16#include <algorithm>
17#include <vector>
18#include <string>
19#include <sstream>
20
21
22#include "BPMaximalMatching.h"
23#include "BPMaximumMatching.h"
24
25using namespace std;
26using namespace combblas;
27
28
33
34// algorithmic options
36int init;
38bool fewexp;
41string ofname;
42
43
44
45/*
46 Remove isolated vertices and purmute
47 */
49{
50
51 int nprocs, myrank;
52 MPI_Comm comm = A.getcommgrid()->GetWorld();
54 MPI_Comm_rank(comm,&myrank);
55
56
59 FullyDistVec<int64_t, int64_t> nonisoRowV; // id's of non-isolated (connected) Row vertices
60 FullyDistVec<int64_t, int64_t> nonisoColV; // id's of non-isolated (connected) Col vertices
61 FullyDistVec<int64_t, int64_t> nonisov; // id's of non-isolated (connected) vertices
62
63 A.Reduce(*ColSums, Column, plus<int64_t>(), static_cast<int64_t>(0));
64 A.Reduce(*RowSums, Row, plus<int64_t>(), static_cast<int64_t>(0));
65
66 // this steps for general graph
67 /*
68 ColSums->EWiseApply(*RowSums, plus<int64_t>()); not needed for bipartite graph
69 nonisov = ColSums->FindInds(bind2nd(greater<int64_t>(), 0));
70 nonisov.RandPerm(); // so that A(v,v) is load-balanced (both memory and time wise)
71 A.operator()(nonisov, nonisov, true); // in-place permute to save memory
72 */
73
74 // this steps for bipartite graph
75 nonisoColV = ColSums->FindInds(bind2nd(greater<int64_t>(), 0));
76 nonisoRowV = RowSums->FindInds(bind2nd(greater<int64_t>(), 0));
77 delete ColSums;
78 delete RowSums;
79
80
81 {
82 nonisoColV.RandPerm();
83 nonisoRowV.RandPerm();
84 }
85
86
87 int64_t nrows1=A.getnrow(), ncols1=A.getncol(), nnz1 = A.getnnz();
88 double avgDeg1 = (double) nnz1/(nrows1+ncols1);
89
90
91 A.operator()(nonisoRowV, nonisoColV, true);
92
93 int64_t nrows2=A.getnrow(), ncols2=A.getncol(), nnz2 = A.getnnz();
94 double avgDeg2 = (double) nnz2/(nrows2+ncols2);
95
96
97 if(myrank == 0)
98 {
99 cout << "ncol nrows nedges deg \n";
100 cout << nrows1 << " " << ncols1 << " " << nnz1 << " " << avgDeg1 << " \n";
101 cout << nrows2 << " " << ncols2 << " " << nnz2 << " " << avgDeg2 << " \n";
102 }
103
105
106
107}
108
109
111{
112 int myrank;
114 if(myrank == 0)
115 {
116 cout << "\n-------------- usage --------------\n";
117 cout << "Usage (random matrix): ./bpmm <er|g500|ssca> <Scale> <EDGEFACTOR> <init><diropt><prune><graft>\n";
118 cout << "Usage (input matrix): ./bpmm <input> <matrix> <init><diropt><prune><graft>\n\n";
119
120 cout << " \n-------------- meaning of arguments ----------\n";
121 cout << "** er: Erdos-Renyi, g500: Graph500 benchmark, ssca: SSCA benchmark\n";
122 cout << "** scale: matrix dimention is 2^scale\n";
123 cout << "** edgefactor: average degree of vertices\n";
124 cout << "** (optional) init : maximal matching algorithm used to initialize\n ";
125 cout << " none: noinit, greedy: greedy init , ks: Karp-Sipser, dmd: dynamic mindegree\n";
126 cout << " default: none\n";
127 cout << "** (optional) randMaximal: random parent selection in greedy/Karp-Sipser\n" ;
128 //cout << "** (optional) diropt: employ direction-optimized BFS\n" ;
129 cout << "** (optional) prune: discard trees as soon as an augmenting path is found\n" ;
130 //cout << "** (optional) graft: employ tree grafting\n" ;
131 cout << "** (optional) moreSplit: more splitting of Matrix.\n" ;
132 cout << "** (optional) randPerm: Randomly permute the matrix for load balance.\n" ;
133 cout << "** (optional) saveMatching: Save the matching vector in a file (filename: inputfile_matching.txt).\n" ;
134 cout << "(order of optional arguments does not matter)\n";
135
136
137 cout << " \n-------------- examples ----------\n";
138 cout << "Example: mpirun -np 4 ./bpmm g500 18 16" << endl;
139 cout << "Example: mpirun -np 4 ./bpmm g500 18 16 ks diropt graft" << endl;
140 cout << "Example: mpirun -np 4 ./bpmm input cage12.mtx randPerm ks diropt graft\n" << endl;
141 }
142}
143
144void GetOptions(char* argv[], int argc)
145{
146 string allArg="";
147 for(int i=0; i<argc; i++)
148 {
149 allArg += string(argv[i]);
150 }
151
152 if(allArg.find("prune")!=string::npos)
153 prune = true;
154 if(allArg.find("fewexp")!=string::npos)
155 fewexp = true;
156 if(allArg.find("moreSplit")!=string::npos)
157 moreSplit = true;
158 if(allArg.find("saveMatching")!=string::npos)
159 saveMatching=true;
160 if(allArg.find("randMM")!=string::npos)
161 randMM = true;
162 if(allArg.find("randMaximal")!=string::npos)
163 randMaximal = true;
164 if(allArg.find("randPerm")!=string::npos)
165 randPerm = true;
166 if(allArg.find("greedy")!=string::npos)
167 init = GREEDY;
168 else if(allArg.find("ks")!=string::npos)
170 else if(allArg.find("dmd")!=string::npos)
171 init = DMD;
172 else
173 init = NO_INIT;
174
175}
176
178{
180 tinfo.str("");
181 tinfo << "\n---------------------------------\n";
182 tinfo << "Calling maximum-cardinality matching with options: " << endl;
183 tinfo << " init: ";
184 if(init == NO_INIT) tinfo << " no-init ";
185 if(init == KARP_SIPSER) tinfo << " Karp-Sipser, ";
186 if(init == DMD) tinfo << " dynamic mindegree, ";
187 if(init == GREEDY) tinfo << " greedy, ";
188 if(randMaximal) tinfo << " random parent selection in greedy/Karp-Sipser, ";
189 if(prune) tinfo << " tree pruning, ";
190 if(moreSplit) tinfo << " moreSplit ";
191 if(randPerm) tinfo << " Randomly permute the matrix for load balance ";
192 if(saveMatching) tinfo << " Write the matcing in a file";
193 tinfo << "\n---------------------------------\n\n";
195
196}
197
199{
200 FullyDistVec<int64_t, int64_t> mateRow2Col ( A.getcommgrid(), A.getnrow(), (int64_t) -1);
201 FullyDistVec<int64_t, int64_t> mateCol2Row ( A.getcommgrid(), A.getncol(), (int64_t) -1);
202
203 // best option
204 init = DMD; randMaximal = false; randMM = true; prune = true;
208 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
209 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
210
211 // best option + KS
212 init = KARP_SIPSER; randMaximal = true; randMM = true; prune = true;
216 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
217 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
218
219
220 // best option + Greedy
221 init = GREEDY; randMaximal = true; randMM = true; prune = true;
225 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
226 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
227
228 // best option + No init
229 init = NO_INIT; randMaximal = false; randMM = true; prune = true;
232 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
233 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
234
235
236 // best option - randMM
237 init = DMD; randMaximal = false; randMM = false; prune = true;
241 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
242 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
243
244
245 // best option - prune
246 init = DMD; randMaximal = false; randMM = true; prune = false;
250 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
251 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
252
253}
254
255
256// default experiment
258{
259 FullyDistVec<int64_t, int64_t> mateRow2Col ( A.getcommgrid(), A.getnrow(), (int64_t) -1);
260 FullyDistVec<int64_t, int64_t> mateCol2Row ( A.getcommgrid(), A.getncol(), (int64_t) -1);
261
262 // best option
263 init = DMD; randMaximal = false; randMM = true; prune = true;
267 if(saveMatching && ofname!="")
268 {
269 mateRow2Col.ParallelWrite(ofname,false,false);
270 }
271 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
272 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
273
274}
275
276
277
278
279
280
281
283{
284
285 int nprocs, myrank;
286 MPI_Comm comm = A.getcommgrid()->GetWorld();
288 MPI_Comm_rank(comm,&myrank);
289
290 FullyDistVec<int64_t, int64_t> mateRow2Col ( A.getcommgrid(), A.getnrow(), (int64_t) -1);
291 FullyDistVec<int64_t, int64_t> mateCol2Row ( A.getcommgrid(), A.getncol(), (int64_t) -1);
292
293 double time_start;
294
295 // best option
296 init = DMD; randMaximal = false; randMM = true; prune = true;
299 double time_dmd = MPI_Wtime()-time_start;
300 int64_t cardDMD = mateRow2Col.Count([](int64_t mate){return mate!=-1;});
301
305 int64_t mmcardDMD = mateRow2Col.Count([](int64_t mate){return mate!=-1;});
306 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
307 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
308
309 // best option + KS
310 init = KARP_SIPSER; randMaximal = true; randMM = true; prune = true;
313 double time_ks = MPI_Wtime()-time_start;
314 int64_t cardKS = mateRow2Col.Count([](int64_t mate){return mate!=-1;});
315
319 int64_t mmcardKS = mateRow2Col.Count([](int64_t mate){return mate!=-1;});
320 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
321 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
322
323
324 // best option + Greedy
325 init = GREEDY; randMaximal = true; randMM = true; prune = true;
329 int64_t cardGreedy = mateRow2Col.Count([](int64_t mate){return mate!=-1;});
330
334 int64_t mmcardGreedy = mateRow2Col.Count([](int64_t mate){return mate!=-1;});
335 mateRow2Col.Apply([](int64_t val){return (int64_t) -1;});
336 mateCol2Row.Apply([](int64_t val){return (int64_t) -1;});
337
338 if(myrank == 0)
339 {
340 cout << "\n maximal matching experiment \n";
341 cout << cardGreedy << " " << mmcardGreedy << " " << time_greedy << " " << time_mm_greedy << " " << cardKS << " " << mmcardKS << " " << time_ks << " " << time_mm_ks << " " << cardDMD << " " << mmcardDMD << " " << time_dmd << " " << time_mm_dmd << " \n";
342 }
343}
344
345
346
347int main(int argc, char* argv[])
348{
349
350 // ------------ initialize MPI ---------------
351 int provided;
354 {
355 printf("ERROR: The MPI library does not have MPI_THREAD_SERIALIZED support\n");
357 }
358 int nprocs, myrank;
361 if(argc < 3)
362 {
363 ShowUsage();
364 MPI_Finalize();
365 return -1;
366 }
367
368 init = DMD;
369 randMaximal = false;
370 prune = false;
371 randMM = true;
372 moreSplit = false;
373 fewexp=false;
374 saveMatching = false;
375 ofname = "";
376 randPerm = false;
377
378 SpParHelper::Print("***** I/O and other preprocessing steps *****\n");
379 // ------------ Process input arguments and build matrix ---------------
380 {
381
384 double t01, t02;
385 if(string(argv[1]) == string("input")) // input option
386 {
387 ABool = new Par_DCSC_Bool();
388
389 string filename(argv[2]);
390 tinfo.str("");
391 tinfo << "\n**** Reading input matrix: " << filename << " ******* " << endl;
393 t01 = MPI_Wtime();
394 ABool->ParallelReadMM(filename, true, maximum<bool>()); // one-based matrix market file
395 t02 = MPI_Wtime();
396 ABool->PrintInfo();
397 tinfo.str("");
398 tinfo << "Reader took " << t02-t01 << " seconds" << endl;
400 GetOptions(argv+3, argc-3);
401 if(saveMatching)
402 {
403 ofname = filename + ".matching.out";
404 }
405
406
407 }
408 else if(argc < 4)
409 {
410 ShowUsage();
411 MPI_Finalize();
412 return -1;
413 }
414 else
415 {
416 unsigned scale = (unsigned) atoi(argv[2]);
417 unsigned EDGEFACTOR = (unsigned) atoi(argv[3]);
418 double initiator[4];
419 if(string(argv[1]) == string("er"))
420 {
421 initiator[0] = .25;
422 initiator[1] = .25;
423 initiator[2] = .25;
424 initiator[3] = .25;
425 if(myrank==0)
426 cout << "Randomly generated ER matric\n";
427 }
428 else if(string(argv[1]) == string("g500"))
429 {
430 initiator[0] = .57;
431 initiator[1] = .19;
432 initiator[2] = .19;
433 initiator[3] = .05;
434 if(myrank==0)
435 cout << "Randomly generated G500 matric\n";
436 }
437 else if(string(argv[1]) == string("ssca"))
438 {
439 initiator[0] = .6;
440 initiator[1] = .4/3;
441 initiator[2] = .4/3;
442 initiator[3] = .4/3;
443 if(myrank==0)
444 cout << "Randomly generated SSCA matric\n";
445 }
446 else
447 {
448 if(myrank == 0)
449 printf("The input type - %s - is not recognized.\n", argv[2]);
451 }
452
453 SpParHelper::Print("Generating input matrix....\n");
454 t01 = MPI_Wtime();
457 ABool = new Par_DCSC_Bool(*DEL, false);
458 delete DEL;
459 t02 = MPI_Wtime();
460 ABool->PrintInfo();
461 tinfo.str("");
462 tinfo << "Generator took " << t02-t01 << " seconds" << endl;
464
466 //removeIsolated(*ABool);
467 SpParHelper::Print("Generated matrix symmetricized....\n");
468 ABool->PrintInfo();
469
470 GetOptions(argv+4, argc-4);
471
472 }
473
474
475 if(randPerm)
476 {
477 // randomly permute for load balance
478 SpParHelper::Print("Performing random permutation of matrix.\n");
481 prow.iota(ABool->getnrow(), 0);
482 pcol.iota(ABool->getncol(), 0);
483 prow.RandPerm();
484 pcol.RandPerm();
485 (*ABool)(prow, pcol, true);
486 SpParHelper::Print("Performed random permutation of matrix.\n");
487 }
488
489
492 AT.Transpose();
493
494 // Reduce is not multithreaded, so I am doing it here
496 A.Reduce(degCol, Column, plus<int64_t>(), static_cast<int64_t>(0));
497
498 // TODO: Follow the AWPM guideline to use CSC matrices.
499 // Currently this file does not use multithreading in SpMSpV
500/*
501 int nthreads;
502#ifdef THREADED
503#pragma omp parallel
504 {
505 int splitPerThread = 1;
506 if(moreSplit) splitPerThread = 4;
507 nthreads = omp_get_num_threads();
508 cblas_splits = nthreads*splitPerThread;
509 }
510 tinfo.str("");
511 tinfo << "Threading activated with " << nthreads << " threads, and matrix split into "<< cblas_splits << " parts" << endl;
512 SpParHelper::Print(tinfo.str());
513 //A.ActivateThreading(cblas_splits); // note: crash on empty matrix
514 //AT.ActivateThreading(cblas_splits);
515#endif
516 */
517
518
519 SpParHelper::Print("**************************************************\n\n");
521 }
522 MPI_Finalize();
523 return 0;
524}
525
526
#define GREEDY
#define KARP_SIPSER
#define DMD
#define NO_INIT
bool moreSplit
void experiment_maximal(Par_DCSC_Bool &A, Par_DCSC_Bool &AT, FullyDistVec< int64_t, int64_t > degCol)
bool randMaximal
SpParMat< int64_t, bool, SpDCCols< int64_t, bool > > Par_DCSC_Bool
bool fewexp
bool prune
SpParMat< int64_t, bool, SpCCols< int64_t, bool > > Par_CSC_Bool
int init
void GetOptions(char *argv[], int argc)
SpParMat< int64_t, double, SpDCCols< int64_t, double > > Par_DCSC_Double
void defaultExp(Par_DCSC_Bool &A, Par_DCSC_Bool &AT, FullyDistVec< int64_t, int64_t > degCol)
void ShowUsage()
bool randMM
string ofname
void showCurOptions()
bool saveMatching
bool randPerm
SpParMat< int64_t, int64_t, SpDCCols< int64_t, int64_t > > Par_DCSC_int64_t
void removeIsolated(Par_DCSC_Bool &A)
void experiment(Par_DCSC_Bool &A, Par_DCSC_Bool &AT, FullyDistVec< int64_t, int64_t > degCol)
int cblas_splits
Definition BcastTest.cpp:26
#define EDGEFACTOR
Definition DirOptBFS.cpp:81
int main()
Definition Driver.cpp:12
void GenGraph500Data(double initiator[4], int log_numverts, int edgefactor, bool scramble=false, bool packed=false)
static void Print(const std::string &s)
int nprocs
Definition comms.cpp:55
long int64_t
Definition compat.h:21
@ Column
Definition SpDefs.h:115
void Symmetricize(PARMAT &A)
Definition ReadMatDist.h:29
void maximumMatching(SpParMat< IT, NT, DER > &A, FullyDistVec< IT, IT > &mateRow2Col, FullyDistVec< IT, IT > &mateCol2Row, bool prune=true, bool randMM=false, bool maximizeWeight=false)
void MaximalMatching(Par_DCSC_Bool &A, Par_DCSC_Bool &AT, FullyDistVec< IT, IT > &mateRow2Col, FullyDistVec< IT, IT > &mateCol2Row, FullyDistVec< IT, IT > &degColRecv, int type, bool rand=true)
double A