COMBINATORIAL_BLAS  1.6
RandomParentBFS.cpp
Go to the documentation of this file.
1 #define DETERMINISTIC
2 #include "CombBLAS/CombBLAS.h"
3 #include <mpi.h>
4 #include <sys/time.h>
5 #include <iostream>
6 #include <functional>
7 #include <algorithm>
8 #include <vector>
9 #include <string>
10 #include <sstream>
11 #ifdef THREADED
12  #ifndef _OPENMP
13  #define _OPENMP
14  #endif
15 #endif
16 
17 #ifdef _OPENMP
18  #include <omp.h>
19 #endif
20 
21 
27 #ifdef _OPENMP
28 int cblas_splits = omp_get_max_threads();
29 #else
30 int cblas_splits = 1;
31 #endif
32 
33 #define ITERS 16
34 #define EDGEFACTOR 16
35 using namespace std;
36 using namespace combblas;
37 
38 
39 MTRand GlobalMT(123); // for reproducable result
40 
41 
42 template <typename PARMAT>
43 void Symmetricize(PARMAT & A)
44 {
45  // boolean addition is practically a "logical or"
46  // therefore this doesn't destruct any links
47  PARMAT AT = A;
48  AT.Transpose();
49  A += AT;
50 }
51 
52 
53 
54 struct ParentType
55 {
56 public:
57  ParentType(){parent=-1; p = 0;};
58  ParentType(int64_t x){parent=(x); p = 0;};
59  friend ostream& operator<<(ostream& os, const ParentType & vertex ){os << "Parent=" << vertex.parent << " p=" << vertex.p; return os;};
60  //private:
62  float p;
63 
64 };
65 
66 
67 
68 
69 
70 
71 struct Edge_randomizer : public std::unary_function<std::pair<bool, float>, std::pair<bool, float>>
72 {
73  const std::pair<bool, float> operator()(const std::pair<bool, float> & x) const
74  {
75  float edgeRand = static_cast<float>(rand()); // random range(0,1)
76  return std::pair<bool, float>(x.first, edgeRand);
77  }
78 };
79 
80 
81 
82 static void MPI_randuniq(void * invec, void * inoutvec, int * len, MPI_Datatype *datatype)
83 {
85  int64_t * inveccast = (int64_t *) invec;
86  int64_t * inoutveccast = (int64_t *) inoutvec;
87  for (int i=0; i<*len; i++ )
88  inoutveccast[i] = RR(inveccast[i], inoutveccast[i]);
89 }
90 
91 
93 {
94  static MPI_Op MPI_BFSRAND;
95  typedef int64_t T_promote;
96  static ParentType id(){ return ParentType(); };
97  static bool returnedSAID() { return false; }
98  //static MPI_Op mpi_op() { return MPI_MAX; }; // do we need this?
99 
100  static ParentType add(const ParentType & arg1, const ParentType & arg2)
101  {
102  //cout << arg1 << " ;;; " << arg2 << endl;
103  if(arg1.p < arg2.p) return arg1;
104  else return arg2;
105  }
106 
107  static ParentType multiply(const T_promote & arg1, const ParentType & arg2)
108  {
109  ParentType temp;
110  temp.parent = arg2.parent;
111  temp.p = GlobalMT.rand();
112  return temp;
113  }
114 
115  static void axpy(T_promote a, const ParentType & x, ParentType & y)
116  {
117  y = add(y, multiply(a, x));
118  }
119 };
120 
121 
122 
123 // This one is used for BFS iteration
125 {
127  static T_promote id(){ return -1; };
128  static bool returnedSAID() { return false; }
129  //static MPI_Op mpi_op() { return MPI_MAX; };
130 
131  static T_promote add(const T_promote & arg1, const T_promote & arg2)
132  {
133  cout << arg1 << " a " << arg2 << endl;
134  return std::max(arg1, arg2);
135  }
136 
137  static T_promote multiply(const bool & arg1, const T_promote & arg2)
138  {
139  cout << arg1 << " m " << arg2 << endl;
140  return arg2;
141  }
142  /*
143  static void axpy(bool a, const T_promote & x, T_promote & y)
144  {
145  y = std::max(y, x);
146  }*/
147 };
148 
152 
153 
155 {
156 
157  int nprocs, myrank;
158  MPI_Comm_size(MPI_COMM_WORLD,&nprocs);
159  MPI_Comm_rank(MPI_COMM_WORLD,&myrank);
160 
161 
163 
164  fringe.SetElement(0, ParentType(0));
165  fringe.SetElement(1, ParentType(1));
166  fringe.SetElement(5, ParentType(5));
167  fringe.SetElement(6, ParentType(6));
168  fringe.SetElement(7, ParentType(7));
169 
170  PSpMat_Int64 A = Aeff;
171  //A.PrintInfo();
172  SpMV<SelectRandSRing>(A, fringe, fringe, false);
173  fringe.DebugPrint();
174 }
175 
176 
177 
178 int main(int argc, char* argv[])
179 {
180  int nprocs, myrank;
181  MPI_Init(&argc, &argv);
182  MPI_Comm_size(MPI_COMM_WORLD,&nprocs);
183  MPI_Comm_rank(MPI_COMM_WORLD,&myrank);
184  if(argc < 2)
185  {
186  if(myrank == 0)
187  {
188  cout << "Usage: ./rpbfs <Scale>" << endl;
189  cout << "Example: mpirun -np 4 ./spbfs 20" << endl;
190  }
191  MPI_Finalize();
192  return -1;
193  }
194  {
195 
196  // Declare objects
197  PSpMat_Bool A;
198  FullyDistVec<int64_t, int64_t> nonisov; // id's of non-isolated (connected) vertices
199  unsigned scale;
200 
201  scale = static_cast<unsigned>(atoi(argv[1]));
202  double initiator[4] = {.57, .19, .19, .05};
203 
204 
206  DEL->GenGraph500Data(initiator, scale, EDGEFACTOR, true, true );
207  MPI_Barrier(MPI_COMM_WORLD);
208 
209 
210  PSpMat_Bool * ABool = new PSpMat_Bool(*DEL, false);
211  delete DEL;
212  int64_t removed = ABool->RemoveLoops();
213  ABool->PrintInfo();
214 
217 
218  ABool->Reduce(*ColSums, Column, plus<int64_t>(), static_cast<int64_t>(0));
219  ABool->Reduce(*RowSums, Row, plus<int64_t>(), static_cast<int64_t>(0));
220  ColSums->EWiseApply(*RowSums, plus<int64_t>());
221  delete RowSums;
222  nonisov = ColSums->FindInds(bind2nd(greater<int64_t>(), 0));
223  delete ColSums;
224  nonisov.RandPerm(); // so that A(v,v) is load-balanced (both memory and time wise)
225 #ifndef NOPERMUTE
226  ABool->operator()(nonisov, nonisov, true); // in-place permute to save memory
227 #endif
228 
229 
230  RandomParentBFS(*ABool);
231 
232 
233  }
234  MPI_Finalize();
235  return 0;
236 }
237 
238 
239 
240 
241 
242 
243 
244 
245 
void SetElement(IT indx, NT numx)
Indexing is performed 0-based.
FullyDistVec< IT, NT > Reduce(Dim dim, _BinaryOperation __binary_op, NT id, _UnaryOperation __unary_op) const
Definition: SpParMat.cpp:791
static bool returnedSAID()
double rand()
static ParentType id()
std::shared_ptr< CommGrid > getcommgrid() const
Definition: SpParMat.h:275
double cblas_alltoalltime
double cblas_transvectime
void GenGraph500Data(double initiator[4], int log_numverts, int edgefactor, bool scramble=false, bool packed=false)
static ParentType add(const ParentType &arg1, const ParentType &arg2)
MTRand GlobalMT(123)
FullyDistVec< IT, IT > FindInds(_Predicate pred) const
Return the indices where pred is true.
double cblas_mergeconttime
static T_promote add(const T_promote &arg1, const T_promote &arg2)
static ParentType multiply(const T_promote &arg1, const ParentType &arg2)
ParentType(int64_t x)
void EWiseApply(const FullyDistVec< IT, NT2 > &other, _BinaryOperation __binary_op, _BinaryPredicate _do_op, const bool useExtendedBinOp)
double cblas_allgathertime
void RandomParentBFS(PSpMat_Bool &Aeff)
With 50/50 chances, return a one of the operants.
Definition: Operations.h:185
double cblas_localspmvtime
void Symmetricize(PARMAT &A)
double A
int cblas_splits
static T_promote multiply(const bool &arg1, const T_promote &arg2)
static T_promote id()
long int64_t
Definition: compat.h:21
#define EDGEFACTOR
static MPI_Op MPI_BFSRAND
const std::pair< bool, float > operator()(const std::pair< bool, float > &x) const
static void axpy(T_promote a, const ParentType &x, ParentType &y)
IT getncol() const
Definition: SpParMat.cpp:694
void PrintInfo() const
Definition: SpParMat.cpp:2387
static bool returnedSAID()
Definition: CCGrid.h:4
SpParMat< int64_t, bool, SpDCCols< int64_t, bool > > PSpMat_Bool
SpParMat< int64_t, int64_t, SpDCCols< int64_t, int64_t > > PSpMat_Int64
SpDCCols< IT, NT > * multiply(SpDCCols< IT, NT > &splitA, SpDCCols< IT, NT > &splitB, CCGrid &CMG, bool isBT, bool threaded)
Definition: Multiplier.h:11
std::ostream & operator<<(std::ostream &os, const TwitterEdge &twe)
Definition: TwitterEdge.h:59
int main(int argc, char *argv[])
SpParMat< int64_t, bool, SpDCCols< int32_t, bool > > PSpMat_s32p64