30 #ifndef __STDC_CONSTANT_MACROS 31 #define __STDC_CONSTANT_MACROS 33 #ifndef __STDC_LIMIT_MACROS 34 #define __STDC_LIMIT_MACROS 61 int main(
int argc,
char* argv[])
64 MPI_Init(&argc, &argv);
65 MPI_Comm_size(MPI_COMM_WORLD,&nprocs);
66 MPI_Comm_rank(MPI_COMM_WORLD,&myrank);
75 cout <<
"Usage: ./betwcent <BASEADDRESS> <K4APPROX> <BATCHSIZE> <output file - optional>" << endl;
76 cout <<
"Example: ./betwcent Data/ 15 128" << endl;
77 cout <<
"Input file input.mtx should be under <BASEADDRESS> in matrix market format" << endl;
78 cout <<
"<BATCHSIZE> should be a multiple of sqrt(p)" << endl;
79 cout <<
"Because <BATCHSIZE> is for the overall matrix (similarly, <K4APPROX> is global as well) " << endl;
86 int K4Approx = atoi(argv[2]);
87 int batchSize = atoi(argv[3]);
89 string directory(argv[1]);
90 string ifilename =
"input.mtx";
91 ifilename = directory+
"/"+ifilename;
93 shared_ptr<CommGrid> fullWorld;
94 fullWorld.reset(
new CommGrid(MPI_COMM_WORLD, 0, 0) );
102 int nPasses = (int) pow(2.0, K4Approx);
103 int numBatches = (int) ceil( static_cast<float>(nPasses)/
static_cast<float>(batchSize));
106 int subBatchSize = batchSize / (AT.
getcommgrid())->GetGridCols();
107 int nBatchSize = subBatchSize * (AT.
getcommgrid())->GetGridCols();
108 nPasses = numBatches * nBatchSize;
110 if(batchSize % (AT.
getcommgrid())->GetGridCols() > 0 && myrank == 0)
112 cout <<
"*** Batchsize is not evenly divisible by the grid dimension ***" << endl;
113 cout <<
"*** Processing "<< nPasses <<
" vertices instead"<< endl;
118 tinfo <<
"Batch processing will occur " << numBatches <<
" times, each processing " << nBatchSize <<
" vertices (overall)" << endl;
121 vector<int> candidates;
125 int nlocpass = numBatches * subBatchSize;
126 while(vertices < nlocpass)
130 single.push_back(vrtxid);
131 int locnnz = ((AT.
seq())(empty,single)).getnnz();
133 MPI_Allreduce( &locnnz, &totnnz, 1, MPI_INT, MPI_SUM, (AT.
getcommgrid())->GetColWorld());
137 candidates.push_back(vrtxid);
144 double t1 = MPI_Wtime();
145 vector<int> batch(subBatchSize);
148 for(
int i=0; i< numBatches; ++i)
150 for(
int j=0; j< subBatchSize; ++j)
152 batch[j] = candidates[i*subBatchSize + j];
160 tuple<int, int, int> * mytuples = NULL;
163 mytuples =
new tuple<int, int, int>[subBatchSize];
164 for(
int k =0; k<subBatchSize; ++k)
166 mytuples[k] = make_tuple(batch[k], k, 1);
176 vector < Dist<bool>::MPI_DCCols * > bfs;
179 while( fringe.
getnnz() > 0 )
183 bfs.push_back(level);
185 fringe = PSpGEMM<PTBOOLINT>(AT, fringe);
192 nspInv.
Apply(bind1st(divides<double>(), 1));
199 for(
int j = bfs.size()-1; j > 0; --j)
205 product =
EWiseMult(product, *bfs[j-1],
false);
206 product =
EWiseMult(product, nsp,
false);
210 for(
int j=0; j < bfs.size(); ++j)
218 bc.
Apply(bind2nd(minus<double>(), nPasses));
220 double t2=MPI_Wtime();
221 double TEPS = (nPasses *
static_cast<float>(A.
getnnz())) / (t2-t1);
224 cout<<
"Computation finished"<<endl;
225 fprintf(stdout,
"%.6lf seconds elapsed for %d starting vertices\n", t2-t1, nPasses);
226 fprintf(stdout,
"TEPS score is: %.6lf\n", TEPS);
228 ofstream output(argv[4]);
SpParMat< int, NT, DCCols > MPI_DCCols
std::shared_ptr< CommGrid > getcommgrid() const
Compute the maximum of two values.
void Apply(_UnaryOperation __unary_op)
DER::LocalIT getlocalrows() const
Dcsc< IU, typename promote_trait< NU1, NU2 >::T_promote > EWiseMult(const Dcsc< IU, NU1 > &A, const Dcsc< IU, NU2 > *B, bool exclude)
void SaveGathered(std::ofstream &outfile, int master, HANDLER handler, bool printProcSplits=false)
DER::LocalIT getlocalcols() const
void Create(const std::vector< IT > &essentials)
static void Print(const std::string &s)
void Apply(_UnaryOperation __unary_op)
SpDCCols< int, NT > DCCols
SpParMat< IT, NT, DER > SubsRefCol(const std::vector< IT > &ci) const
Column indexing with special parallel semantics.
int main(int argc, char *argv[])
void ParallelReadMM(const std::string &filename, bool onebased, _BinaryOperation BinOp)
void EWiseScale(const DenseParMat< IT, NT > &rhs)