COMBINATORIAL_BLAS  1.6
MultTest.cpp
Go to the documentation of this file.
1 /****************************************************************/
2 /* Parallel Combinatorial BLAS Library (for Graph Computations) */
3 /* version 1.6 -------------------------------------------------*/
4 /* date: 6/15/2017 ---------------------------------------------*/
5 /* authors: Ariful Azad, Aydin Buluc --------------------------*/
6 /****************************************************************/
7 /*
8  Copyright (c) 2010-2017, The Regents of the University of California
9 
10  Permission is hereby granted, free of charge, to any person obtaining a copy
11  of this software and associated documentation files (the "Software"), to deal
12  in the Software without restriction, including without limitation the rights
13  to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
14  copies of the Software, and to permit persons to whom the Software is
15  furnished to do so, subject to the following conditions:
16 
17  The above copyright notice and this permission notice shall be included in
18  all copies or substantial portions of the Software.
19 
20  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
21  IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
22  FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
23  AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
24  LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
25  OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
26  THE SOFTWARE.
27  */
28 
29 #include <mpi.h>
30 #include <sys/time.h>
31 #include <iostream>
32 #include <functional>
33 #include <algorithm>
34 #include <vector>
35 #include <sstream>
36 #include "CombBLAS/CombBLAS.h"
37 
38 using namespace std;
39 using namespace combblas;
40 
41 #ifdef TIMING
42 double cblas_alltoalltime;
43 double cblas_allgathertime;
44 #endif
45 
46 #ifdef _OPENMP
47 int cblas_splits = omp_get_max_threads();
48 #else
49 int cblas_splits = 1;
50 #endif
51 
52 
53 // Simple helper class for declarations: Just the numerical type is templated
54 // The index type and the sequential matrix type stays the same for the whole code
55 // In this case, they are "int" and "SpDCCols"
56 template <class NT>
57 class PSpMat
58 {
59 public:
62 };
63 
64 int main(int argc, char* argv[])
65 {
66  int nprocs, myrank;
67  MPI_Init(&argc, &argv);
68  MPI_Comm_size(MPI_COMM_WORLD,&nprocs);
69  MPI_Comm_rank(MPI_COMM_WORLD,&myrank);
70 
71  if(argc < 6)
72  {
73  if(myrank == 0)
74  {
75  cout << "Usage: ./MultTest <MatrixA> <MatrixB> <MatrixC> <vecX> <vecY>" << endl;
76  cout << "<MatrixA>,<MatrixB>,<MatrixC> are absolute addresses, and files should be in triples format" << endl;
77  }
78  MPI_Finalize();
79  return -1;
80  }
81  {
82  string Aname(argv[1]);
83  string Bname(argv[2]);
84  string Cname(argv[3]);
85  string V1name(argv[4]);
86  string V2name(argv[5]);
87 
88  ifstream vecinpx(V1name.c_str());
89  ifstream vecinpy(V2name.c_str());
90 
91  MPI_Barrier(MPI_COMM_WORLD);
92  typedef PlusTimesSRing<double, double> PTDOUBLEDOUBLE;
94 
95  shared_ptr<CommGrid> fullWorld;
96  fullWorld.reset( new CommGrid(MPI_COMM_WORLD, 0, 0) );
97 
98  // construct objects
99  PSpMat<double>::MPI_DCCols A(fullWorld);
100  PSpMat<double>::MPI_DCCols B(fullWorld);
101  PSpMat<double>::MPI_DCCols C(fullWorld);
102  PSpMat<double>::MPI_DCCols CControl(fullWorld);
103  FullyDistVec<int64_t, double> ycontrol(fullWorld);
104  FullyDistVec<int64_t, double> x(fullWorld);
105  FullyDistSpVec<int64_t, double> spycontrol(fullWorld);
106  FullyDistSpVec<int64_t, double> spx(fullWorld);
107 
108  A.ParallelReadMM(Aname, true, maximum<double>());
109 #ifndef NOGEMM
110  B.ParallelReadMM(Bname, true, maximum<double>());
111 
112  CControl.ParallelReadMM(Cname, true, maximum<double>());
113 #endif
114  x.ReadDistribute(vecinpx, 0);
115  spx.ReadDistribute(vecinpx, 0);
116  ycontrol.ReadDistribute(vecinpy,0);
117  spycontrol.ReadDistribute(vecinpy,0);
118 
119  FullyDistVec<int64_t, double> y = SpMV<PTDOUBLEDOUBLE>(A, x);
120  if (ycontrol == y)
121  {
122  SpParHelper::Print("Dense SpMV (fully dist) working correctly\n");
123  }
124  else
125  {
126  SpParHelper::Print("ERROR in Dense SpMV, go fix it!\n");
127  y.ParallelWrite("ycontrol_dense.txt",true);
128  }
129 
130  //FullyDistSpVec<int64_t, double> spy = SpMV<PTDOUBLEDOUBLE>(A, spx);
131 
132  FullyDistSpVec<int64_t, double> spy(spx.getcommgrid(), A.getnrow());
133  SpMV<PTDOUBLEDOUBLE>(A, spx, spy, false);
134 
135  if (spycontrol == spy)
136  {
137  SpParHelper::Print("Sparse SpMV (fully dist) working correctly\n");
138  }
139  else
140  {
141  SpParHelper::Print("ERROR in Sparse SpMV, go fix it!\n");
142  spy.ParallelWrite("ycontrol_sparse.txt",true);
143  }
144 
145  // Test SpMSpV-bucket for general CSC matrices
148  FullyDistSpVec<int64_t, double> spy_csc(spx.getcommgrid(), ACsc.getnrow());
149  SpMV<PTDOUBLEDOUBLE>(ACsc, spx, spy_csc, false, SPA);
150 
151  if (spy == spy_csc)
152  {
153  SpParHelper::Print("SpMSpV-bucket works correctly for general CSC matrices\n");
154  }
155  else
156  {
157  SpParHelper::Print("SpMSpV-bucket does not work correctly for general CSC matrices, go fix it!\n");
158  }
159 
160 
161 #ifndef NOGEMM
162  C = Mult_AnXBn_Synch<PTDOUBLEDOUBLE, double, PSpMat<double>::DCCols >(A,B);
163  if (CControl == C)
164  {
165  SpParHelper::Print("Synchronous Multiplication working correctly\n");
166  // C.SaveGathered("CControl.txt");
167  }
168  else
169  {
170  SpParHelper::Print("ERROR in Synchronous Multiplication, go fix it!\n");
171  }
172 
173  C = Mult_AnXBn_DoubleBuff<PTDOUBLEDOUBLE, double, PSpMat<double>::DCCols >(A,B);
174  if (CControl == C)
175  {
176  SpParHelper::Print("Double buffered multiplication working correctly\n");
177  }
178  else
179  {
180  SpParHelper::Print("ERROR in double buffered multiplication, go fix it!\n");
181  }
182 #endif
184  PSpMat<bool>::MPI_DCCols ABool(A);
185 
186  spx.Apply(bind1st (multiplies<double>(), 100));
187  FullyDistSpVec<int64_t, int64_t> spxint64 (spx);
188  //FullyDistSpVec<int64_t, int64_t> spyint64 = SpMV<SR>(ABool, spxint64, false);
189  FullyDistSpVec<int64_t, int64_t> spyint64(spxint64.getcommgrid(), ABool.getnrow());
190  SpMV<SR>(ABool, spxint64, spyint64, false);
191 
192 
193  ABool.OptimizeForGraph500(optbuf);
194  //FullyDistSpVec<int64_t, int64_t> spyint64buf = SpMV<SR>(ABool, spxint64, false, optbuf);
195  FullyDistSpVec<int64_t, int64_t> spyint64buf(spxint64.getcommgrid(), ABool.getnrow());
196  SpMV<SR>(ABool, spxint64, spyint64buf, false, optbuf);
197 
198 
199  if (spyint64 == spyint64buf)
200  {
201  SpParHelper::Print("Graph500 Optimizations are correct\n");
202  }
203  else
204  {
205  SpParHelper::Print("ERROR in graph500 optimizations, go fix it!\n");
206  spyint64.ParallelWrite("Original_SpMSV.txt",true);
207  spyint64buf.ParallelWrite("Buffered_SpMSV.txt",true);
208  }
209 
210 
211  ABool.ActivateThreading(cblas_splits);
212  //FullyDistSpVec<int64_t, int64_t> spyint64_threaded = SpMV<SR>(ABool, spxint64, false);
213  FullyDistSpVec<int64_t, int64_t> spyint64_threaded(spxint64.getcommgrid(), ABool.getnrow());
214  SpMV<SR>(ABool, spxint64, spyint64_threaded, false);
215 
216  if (spyint64 == spyint64_threaded)
217  {
218  SpParHelper::Print("Multithreaded Sparse SpMV works\n");
219  }
220  else
221  {
222  SpParHelper::Print("ERROR in multithreaded sparse SpMV, go fix it!\n");
223  }
224 
225 
226  // Test SpMSpV-bucket for Boolean CSC matrices
228  PreAllocatedSPA<int64_t> SPA1(ABoolCsc.seq(), cblas_splits*4);
229  FullyDistSpVec<int64_t, int64_t> spyint64_csc_threaded(spxint64.getcommgrid(), ABoolCsc.getnrow());
230  SpMV<SR>(ABoolCsc, spxint64, spyint64_csc_threaded, false, SPA1);
231 
232  if (spyint64 == spyint64_csc_threaded)
233  {
234  SpParHelper::Print("SpMSpV-bucket works correctly for Boolean CSC matrices\n");
235  }
236  else
237  {
238  SpParHelper::Print("ERROR in SpMSpV-bucket with Boolean CSC matrices, go fix it!\n");
239  }
240 
241 
242  vecinpx.clear();
243  vecinpx.close();
244  vecinpy.clear();
245  vecinpy.close();
246  }
247  MPI_Finalize();
248  return 0;
249 }
250 
double B
SpDCCols< int64_t, NT > DCCols
Definition: MultTest.cpp:60
Compute the maximum of two values.
Definition: Operations.h:154
int cblas_splits
Definition: MultTest.cpp:49
SelectMaxSRing< bool, int32_t > SR
Definition: SpMSpVBench.cpp:59
std::shared_ptr< CommGrid > getcommgrid() const
double A
double cblas_allgathertime
Definition: DirOptBFS.cpp:58
std::ifstream & ReadDistribute(std::ifstream &infile, int master, HANDLER handler)
Totally obsolete version that only accepts an ifstream object and ascii files.
void ParallelWrite(const std::string &filename, bool onebased, HANDLER handler, bool includeindices=true)
Definition: FullyDistVec.h:96
double C
void ParallelWrite(const std::string &filename, bool onebased, HANDLER handler, bool includeindices=true, bool includeheader=false)
int main(int argc, char *argv[])
Definition: MultTest.cpp:64
SpParMat< int64_t, NT, DCCols > MPI_DCCols
Definition: MultTest.cpp:61
Definition: CCGrid.h:4
double cblas_alltoalltime
Definition: DirOptBFS.cpp:57
std::ifstream & ReadDistribute(std::ifstream &infile, int master, HANDLER handler)
void ParallelReadMM(const std::string &filename, bool onebased, _BinaryOperation BinOp)
Definition: SpParMat.cpp:3417
IT getnrow() const
Definition: SpParMat.cpp:685