COMBINATORIAL_BLAS  1.6
IndexingTiming.cpp
Go to the documentation of this file.
1 #include <mpi.h>
2 #include <sys/time.h>
3 #include <iostream>
4 #include <functional>
5 #include <algorithm>
6 #include <vector>
7 #include <sstream>
8 #include "CombBLAS/CombBLAS.h"
9 
10 using namespace std;
11 using namespace combblas;
12 
13 #define ITERATIONS 10
14 
15 template <class T>
16 bool from_string(T & t, const string& s, std::ios_base& (*f)(std::ios_base&))
17 {
18  istringstream iss(s);
19  return !(iss >> f >> t).fail();
20 }
21 
22 int main(int argc, char* argv[])
23 {
24  int nprocs, myrank;
25  MPI_Init(&argc, &argv);
26  MPI_Comm_size(MPI_COMM_WORLD,&nprocs);
27  MPI_Comm_rank(MPI_COMM_WORLD,&myrank);
28 
29  if(argc < 3)
30  {
31  if(myrank == 0)
32  {
33  cout << "Usage: ./IndexingTiming Input/Force/Binary <Inputfile>/<Scale>/<BinaryFile>" << endl;
34  }
35  MPI_Finalize();
36  return -1;
37  }
38  {
40  PARDBMAT * A; // declare objects
41  if(string(argv[1]) == string("Input"))
42  {
43  A->ReadDistribute(argv[2], 0);
44  }
45  else if(string(argv[1]) == string("Binary"))
46  {
47  uint64_t n, m;
48  from_string(n,string(argv[3]),std::dec);
49  from_string(m,string(argv[4]),std::dec);
50 
51  ostringstream outs;
52  outs << "Reading " << argv[2] << " with " << n << " vertices and " << m << " edges" << endl;
53  SpParHelper::Print(outs.str());
54  DistEdgeList<int64_t> * DEL = new DistEdgeList<int64_t>(argv[2], n, m);
55  SpParHelper::Print("Read binary input to distributed edge list\n");
56 
57  PermEdges(*DEL);
58  SpParHelper::Print("Permuted Edges\n");
59 
60  RenameVertices(*DEL);
61  SpParHelper::Print("Renamed Vertices\n");
62 
63  A = new PARDBMAT(*DEL, false);
64  delete DEL; // free memory before symmetricizing
65  SpParHelper::Print("Created double Sparse Matrix\n");
66  }
67  else if(string(argv[1]) == string("Force"))
68  {
69  double initiator[4] = {.6, .4/3, .4/3, .4/3};
71 
72  int scale = static_cast<unsigned>(atoi(argv[2]));
73  ostringstream outs;
74  outs << "Forcing scale to : " << scale << endl;
75  SpParHelper::Print(outs.str());
76  DEL->GenGraph500Data(initiator, scale, 8 * ((int64_t) std::pow(2.0, (double) scale)) / nprocs );
77  SpParHelper::Print("Generated local RMAT matrices\n");
78 
79  PermEdges(*DEL);
80  SpParHelper::Print("Permuted Edges\n");
81 
82  RenameVertices(*DEL);
83  SpParHelper::Print("Renamed Vertices\n");
84 
85  // conversion from distributed edge list, keeps self-loops, sums duplicates
86  A = new PARDBMAT(*DEL, false);
87  delete DEL; // free memory before symmetricizing
88  SpParHelper::Print("Created double Sparse Matrix\n");
89  }
90 
91 
92  A->PrintInfo();
94  p.iota(A->getnrow(), 0);
95  p.RandPerm();
96  SpParHelper::Print("Permutation Generated\n");
97  PARDBMAT B = (*A)(p,p);
98  B.PrintInfo();
99 
100  float oldbalance = A->LoadImbalance();
101  float newbalance = B.LoadImbalance();
102  ostringstream outs;
103  outs << "Running on " << nprocs << " cores" << endl;
104  outs << "Old balance: " << oldbalance << endl;
105  outs << "New balance: " << newbalance << endl;
106  SpParHelper::Print(outs.str());
107 
108  MPI_Barrier(MPI_COMM_WORLD);
109  double t1 = MPI_Wtime(); // initilize (wall-clock) timer
110 
111  for(int i=0; i<ITERATIONS; i++)
112  {
113  B = (*A)(p,p);
114  }
115 
116  MPI_Barrier(MPI_COMM_WORLD);
117  double t2 = MPI_Wtime();
118 
119  if(myrank == 0)
120  {
121  cout<<"Indexing Iterations finished"<<endl;
122  printf("%.6lf seconds elapsed per iteration\n", (t2-t1)/(double)ITERATIONS);
123  }
124 
125  // Test #2
126  int nclust = 10;
127  vector< FullyDistVec<int,int> > clusters(nclust);
128  int nperclus = A->getnrow() / nclust;
129 
130  for(int i = 0; i< nclust; i++)
131  {
132  int k = std::min(nperclus, A->getnrow() - nperclus * i);
133  clusters[i].iota(k, nperclus * i);
134  clusters[i] = p(clusters[i]);
135  }
136 
137  for(int i=0; i< nclust; i++)
138  {
139  B = (*A)(clusters[i], clusters[i]);
140  B.PrintInfo();
141  }
142 
143  MPI_Barrier(MPI_COMM_WORLD);
144  t1 = MPI_Wtime(); // initilize (wall-clock) timer
145  for(int i=0; i< nclust; i++)
146  {
147  for(int j=0; j < ITERATIONS; j++)
148  B = (*A)(clusters[i], clusters[i]);
149  }
150 
151  MPI_Barrier(MPI_COMM_WORLD);
152  t2 = MPI_Wtime();
153 
154  if(myrank == 0)
155  {
156  cout<<"Indexing Iterations finished"<<endl;
157  printf("%.6lf seconds elapsed per iteration\n", (t2-t1)/(double)ITERATIONS);
158  }
159 
160  // Test #3: Pruning for SpAsgn
161  for(int i=0; i< nclust; i++)
162  {
163  PARDBMAT C = *A;
164  C.Prune(clusters[i], clusters[i]);
165  C.PrintInfo();
166 
167  B = (*A)(clusters[i], clusters[i]);
168  C.SpAsgn(clusters[i], clusters[i], B);
169 
170  // We should get the original A back.
171  if( C == (*A))
172  {
173  SpParHelper::Print("Pruning and SpAsgn seem to be working\n");
174  }
175  else
176  {
177  SpParHelper::Print("A and C don't match, below is C's info followed by B's info\n");
178  C.PrintInfo();
179  B.PrintInfo();
180  }
181  }
182 
183 
184  delete A;
185  }
186  MPI_Finalize();
187  return 0;
188 }
double B
void GenGraph500Data(double initiator[4], int log_numverts, int edgefactor, bool scramble=false, bool packed=false)
#define ITERATIONS
void PermEdges(DistEdgeList< IT > &DEL)
void RenameVertices(DistEdgeList< IU > &DEL)
double A
void iota(IT globalsize, NT first)
long int64_t
Definition: compat.h:21
double C
int main(int argc, char *argv[])
Definition: CCGrid.h:4
bool from_string(T &t, const string &s, std::ios_base &(*f)(std::ios_base &))