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 30 #ifndef _SP_DCCOLS_H 31 #define _SP_DCCOLS_H 32 33 #include <iostream> 34 #include <fstream> 35 #include <cmath> 36 #include "SpMat.h" // Best to include the base class first 37 #include "SpHelper.h" 38 #include "StackEntry.h" 39 #include "dcsc.h" 40 #include "Isect.h" 41 #include "Semirings.h" 42 #include "MemoryPool.h" 43 #include "LocArr.h" 44 #include "Friends.h" 45 #include "CombBLAS.h" 46 #include "FullyDist.h" 47 48 namespace combblas { 49 50 template <class IT, class NT> 51 class SpDCCols: public SpMat<IT, NT, SpDCCols<IT, NT> > 52 { 53 public: 54 typedef IT LocalIT; 55 typedef NT LocalNT; 56 57 // Constructors : 58 SpDCCols (); 59 SpDCCols (IT size, IT nRow, IT nCol, IT nzc); 60 SpDCCols (const SpTuples<IT,NT> & rhs, bool transpose); 61 SpDCCols (IT nRow, IT nCol, IT nnz1, const std::tuple<IT, IT, NT> * rhs, bool transpose); 62 63 SpDCCols (const SpDCCols<IT,NT> & rhs); // Actual copy constructor 64 ~SpDCCols(); 65 66 template <typename NNT> operator SpDCCols<IT,NNT> () const; //!< NNT: New numeric type 67 template <typename NIT, typename NNT> operator SpDCCols<NIT,NNT> () const; //!< NNT: New numeric type, NIT: New index type 68 69 // Member Functions and Operators: 70 SpDCCols<IT,NT> & operator= (const SpDCCols<IT, NT> & rhs); 71 SpDCCols<IT,NT> & operator+= (const SpDCCols<IT, NT> & rhs); 72 SpDCCols<IT,NT> operator() (IT ri, IT ci) const; 73 SpDCCols<IT,NT> operator() (const std::vector<IT> & ri, const std::vector<IT> & ci) const; 74 bool operator== (const SpDCCols<IT, NT> & rhs) const 75 { 76 if(rhs.nnz == 0 && nnz == 0) 77 return true; 78 if(nnz != rhs.nnz || m != rhs.m || n != rhs.n) 79 return false; 80 return ((*dcsc) == (*(rhs.dcsc))); 81 } 82 83 84 class SpColIter //! Iterate over (sparse) columns of the sparse matrix 85 { 86 public: 87 class NzIter //! Iterate over the nonzeros of the sparse column 88 { 89 public: rid(ir)90 NzIter(IT * ir = NULL, NT * num = NULL) : rid(ir), val(num) {} 91 92 bool operator==(const NzIter & other) 93 { 94 return(rid == other.rid); // compare pointers 95 } 96 bool operator!=(const NzIter & other) 97 { 98 return(rid != other.rid); 99 } 100 bool operator<(const NzIter & other) 101 { 102 return(rid < other.rid); 103 } 104 NzIter operator+(IT inc) 105 { 106 rid+=inc; 107 val+=inc; 108 return(*this); 109 } 110 NzIter operator-(IT inc) 111 { 112 rid-=inc; 113 val-=inc; 114 return(*this); 115 } 116 NzIter & operator++() // prefix operator 117 { 118 ++rid; 119 ++val; 120 return(*this); 121 } 122 NzIter operator++(int) // postfix operator 123 { 124 NzIter tmp(*this); 125 ++(*this); 126 return(tmp); 127 } rowid()128 IT rowid() const //!< Return the "local" rowid of the current nonzero entry. 129 { 130 return (*rid); 131 } value()132 NT & value() //!< value is returned by reference for possible updates 133 { 134 return (*val); 135 } 136 private: 137 IT * rid; 138 NT * val; 139 140 }; 141 cptr(cp)142 SpColIter(IT * cp = NULL, IT * jc = NULL) : cptr(cp), cid(jc) {} 143 bool operator==(const SpColIter& other) 144 { 145 return(cptr == other.cptr); // compare pointers 146 } 147 bool operator!=(const SpColIter& other) 148 { 149 return(cptr != other.cptr); 150 } 151 152 SpColIter& operator++() // prefix operator 153 { 154 ++cptr; 155 ++cid; 156 return(*this); 157 } 158 SpColIter operator++(int) // postfix operator (common) 159 { 160 SpColIter tmp(*this); 161 ++(*this); 162 return(tmp); 163 } 164 SpColIter operator+(IT inc) 165 { 166 cptr+=inc; 167 cid+=inc; 168 return(*this); 169 } 170 SpColIter operator-(IT inc) 171 { 172 cptr-=inc; 173 cid-=inc; 174 return(*this); 175 } colid()176 IT colid() const //!< Return the "local" colid of the current column. 177 { 178 return (*cid); 179 } colptr()180 IT colptr() const 181 { 182 return (*cptr); 183 } colptrnext()184 IT colptrnext() const 185 { 186 return (*(cptr+1)); 187 } nnz()188 IT nnz() const 189 { 190 return (colptrnext() - colptr()); 191 } 192 private: 193 IT * cptr; 194 IT * cid; 195 }; 196 begcol()197 SpColIter begcol() 198 { 199 if( nnz > 0 ) 200 return SpColIter(dcsc->cp, dcsc->jc); 201 else 202 return SpColIter(NULL, NULL); 203 } begcol(int i)204 SpColIter begcol(int i) // multithreaded version 205 { 206 if( dcscarr[i] ) 207 return SpColIter(dcscarr[i]->cp, dcscarr[i]->jc); 208 else 209 return SpColIter(NULL, NULL); 210 } 211 endcol()212 SpColIter endcol() 213 { 214 if( nnz > 0 ) 215 return SpColIter(dcsc->cp + dcsc->nzc, NULL); 216 else 217 return SpColIter(NULL, NULL); 218 } 219 endcol(int i)220 SpColIter endcol(int i) //multithreaded version 221 { 222 if( dcscarr[i] ) 223 return SpColIter(dcscarr[i]->cp + dcscarr[i]->nzc, NULL); 224 else 225 return SpColIter(NULL, NULL); 226 } 227 begnz(const SpColIter & ccol)228 typename SpColIter::NzIter begnz(const SpColIter & ccol) //!< Return the beginning iterator for the nonzeros of the current column 229 { 230 return typename SpColIter::NzIter( dcsc->ir + ccol.colptr(), dcsc->numx + ccol.colptr() ); 231 } 232 endnz(const SpColIter & ccol)233 typename SpColIter::NzIter endnz(const SpColIter & ccol) //!< Return the ending iterator for the nonzeros of the current column 234 { 235 return typename SpColIter::NzIter( dcsc->ir + ccol.colptrnext(), NULL ); 236 } 237 begnz(const SpColIter & ccol,int i)238 typename SpColIter::NzIter begnz(const SpColIter & ccol, int i) //!< multithreaded version 239 { 240 return typename SpColIter::NzIter( dcscarr[i]->ir + ccol.colptr(), dcscarr[i]->numx + ccol.colptr() ); 241 } 242 endnz(const SpColIter & ccol,int i)243 typename SpColIter::NzIter endnz(const SpColIter & ccol, int i) //!< multithreaded version 244 { 245 return typename SpColIter::NzIter( dcscarr[i]->ir + ccol.colptrnext(), NULL ); 246 } 247 248 template <typename _UnaryOperation> Apply(_UnaryOperation __unary_op)249 void Apply(_UnaryOperation __unary_op) 250 { 251 if(nnz > 0) 252 dcsc->Apply(__unary_op); 253 } 254 255 template <typename _UnaryOperation, typename GlobalIT> 256 SpDCCols<IT,NT>* PruneI(_UnaryOperation __unary_op, bool inPlace, GlobalIT rowOffset, GlobalIT colOffset); 257 template <typename _UnaryOperation> 258 SpDCCols<IT,NT>* Prune(_UnaryOperation __unary_op, bool inPlace); 259 template <typename _BinaryOperation> 260 SpDCCols<IT,NT>* PruneColumn(NT* pvals, _BinaryOperation __binary_op, bool inPlace); 261 template <typename _BinaryOperation> 262 SpDCCols<IT,NT>* PruneColumn(IT* pinds, NT* pvals, _BinaryOperation __binary_op, bool inPlace); 263 264 template <typename _BinaryOperation> UpdateDense(NT ** array,_BinaryOperation __binary_op)265 void UpdateDense(NT ** array, _BinaryOperation __binary_op) const 266 { 267 if(nnz > 0 && dcsc != NULL) 268 dcsc->UpdateDense(array, __binary_op); 269 } 270 271 void EWiseScale(NT ** scaler, IT m_scaler, IT n_scaler); 272 void EWiseMult (const SpDCCols<IT,NT> & rhs, bool exclude); 273 274 void Transpose(); //!< Mutator version, replaces the calling object 275 SpDCCols<IT,NT> TransposeConst() const; //!< Const version, doesn't touch the existing object 276 SpDCCols<IT,NT> * TransposeConstPtr() const; 277 RowSplit(int numsplits)278 void RowSplit(int numsplits) 279 { 280 BooleanRowSplit(*this, numsplits); // only works with boolean arrays 281 } 282 283 void ColSplit(int parts, std::vector< SpDCCols<IT,NT> > & matrices); //!< \attention Destroys calling object (*this) 284 void ColConcatenate(std::vector< SpDCCols<IT,NT> > & matrices); //!< \attention Destroys its parameters (matrices) 285 286 void Split(SpDCCols<IT,NT> & partA, SpDCCols<IT,NT> & partB); //!< \attention Destroys calling object (*this) 287 void Merge(SpDCCols<IT,NT> & partA, SpDCCols<IT,NT> & partB); //!< \attention Destroys its parameters (partA & partB) 288 289 void CreateImpl(const std::vector<IT> & essentials); 290 void CreateImpl(IT size, IT nRow, IT nCol, std::tuple<IT, IT, NT> * mytuples); 291 void CreateImpl(IT * _cp, IT * _jc, IT * _ir, NT * _numx, IT _nz, IT _nzc, IT _m, IT _n); 292 293 294 Arr<IT,NT> GetArrays() const; 295 std::vector<IT> GetEssentials() const; 296 const static IT esscount; 297 isZero()298 bool isZero() const { return (nnz == 0); } getnrow()299 IT getnrow() const { return m; } getncol()300 IT getncol() const { return n; } getnnz()301 IT getnnz() const { return nnz; } getnzc()302 IT getnzc() const { return (nnz == 0) ? 0: dcsc->nzc; } getnsplit()303 int getnsplit() const { return splits; } 304 305 std::ofstream& put(std::ofstream & outfile) const; 306 std::ifstream& get(std::ifstream & infile); 307 void PrintInfo() const; 308 void PrintInfo(std::ofstream & out) const; 309 310 template <typename SR> 311 int PlusEq_AtXBt(const SpDCCols<IT,NT> & A, const SpDCCols<IT,NT> & B); 312 313 template <typename SR> 314 int PlusEq_AtXBn(const SpDCCols<IT,NT> & A, const SpDCCols<IT,NT> & B); 315 316 template <typename SR> 317 int PlusEq_AnXBt(const SpDCCols<IT,NT> & A, const SpDCCols<IT,NT> & B); 318 319 template <typename SR> 320 int PlusEq_AnXBn(const SpDCCols<IT,NT> & A, const SpDCCols<IT,NT> & B); 321 GetDCSC()322 Dcsc<IT, NT> * GetDCSC() const // only for single threaded matrices 323 { 324 return dcsc; 325 } 326 GetDCSC(int i)327 Dcsc<IT, NT> * GetDCSC(int i) const // only for split (multithreaded) matrices 328 { 329 return dcscarr[i]; 330 } 331 GetInternal()332 auto GetInternal() const { return GetDCSC(); } GetInternal(int i)333 auto GetInternal(int i) const { return GetDCSC(i); } 334 335 336 private: 337 void CopyDcsc(Dcsc<IT,NT> * source); 338 SpDCCols<IT,NT> ColIndex(const std::vector<IT> & ci) const; //!< col indexing without multiplication 339 340 template <typename SR, typename NTR> 341 SpDCCols< IT, typename promote_trait<NT,NTR>::T_promote > OrdOutProdMult(const SpDCCols<IT,NTR> & rhs) const; 342 343 template <typename SR, typename NTR> 344 SpDCCols< IT, typename promote_trait<NT,NTR>::T_promote > OrdColByCol(const SpDCCols<IT,NTR> & rhs) const; 345 346 SpDCCols (IT size, IT nRow, IT nCol, const std::vector<IT> & indices, bool isRow); // Constructor for indexing 347 SpDCCols (IT nRow, IT nCol, Dcsc<IT,NT> * mydcsc); // Constructor for multiplication 348 349 // Anonymous union 350 union { 351 Dcsc<IT, NT> * dcsc; 352 Dcsc<IT, NT> ** dcscarr; 353 }; 354 355 IT m; 356 IT n; 357 IT nnz; 358 359 int splits; // for multithreading 360 361 template <class IU, class NU> 362 friend class SpDCCols; // Let other template instantiations (of the same class) access private members 363 364 template <class IU, class NU> 365 friend class SpTuples; 366 367 // AL: removed this because it appears illegal and causes this compiler warning: 368 // warning: dependent nested name specifier 'SpDCCols<IU, NU>::' for friend class declaration is not supported; turning off access control for 'SpDCCols' 369 //template <class IU, class NU> 370 //friend class SpDCCols<IU, NU>::SpColIter; 371 372 template<typename IU> 373 friend void BooleanRowSplit(SpDCCols<IU, bool> & A, int numsplits); 374 375 template<typename IU, typename NU1, typename NU2> 376 friend SpDCCols<IU, typename promote_trait<NU1,NU2>::T_promote > EWiseMult (const SpDCCols<IU,NU1> & A, const SpDCCols<IU,NU2> & B, bool exclude); 377 378 template<typename N_promote, typename IU, typename NU1, typename NU2, typename _BinaryOperation> 379 friend SpDCCols<IU, N_promote > EWiseApply (const SpDCCols<IU,NU1> & A, const SpDCCols<IU,NU2> & B, _BinaryOperation __binary_op, bool notB, const NU2& defaultBVal); 380 381 template <typename RETT, typename IU, typename NU1, typename NU2, typename _BinaryOperation, typename _BinaryPredicate> 382 friend SpDCCols<IU,RETT> EWiseApply (const SpDCCols<IU,NU1> & A, const SpDCCols<IU,NU2> & B, _BinaryOperation __binary_op, _BinaryPredicate do_op, bool allowANulls, bool allowBNulls, const NU1& ANullVal, const NU2& BNullVal, const bool allowIntersect); 383 384 template<class SR, class NUO, class IU, class NU1, class NU2> 385 friend SpTuples<IU, NUO> * Tuples_AnXBn 386 (const SpDCCols<IU, NU1> & A, const SpDCCols<IU, NU2> & B, bool clearA, bool clearB); 387 388 template<class SR, class NUO, class IU, class NU1, class NU2> 389 friend SpTuples<IU, NUO> * Tuples_AnXBt 390 (const SpDCCols<IU, NU1> & A, const SpDCCols<IU, NU2> & B, bool clearA, bool clearB); 391 392 template<class SR, class NUO, class IU, class NU1, class NU2> 393 friend SpTuples<IU, NUO> * Tuples_AtXBn 394 (const SpDCCols<IU, NU1> & A, const SpDCCols<IU, NU2> & B, bool clearA, bool clearB); 395 396 template<class SR, class NUO, class IU, class NU1, class NU2> 397 friend SpTuples<IU, NUO> * Tuples_AtXBt 398 (const SpDCCols<IU, NU1> & A, const SpDCCols<IU, NU2> & B, bool clearA, bool clearB); 399 400 template <typename SR, typename IU, typename NU, typename RHS, typename LHS> 401 friend void dcsc_gespmv (const SpDCCols<IU, NU> & A, const RHS * x, LHS * y); 402 403 template <typename SR, typename IU, typename NU, typename RHS, typename LHS> 404 friend void dcsc_gespmv_threaded (const SpDCCols<IU, NU> & A, const RHS * x, LHS * y); 405 406 template <typename SR, typename IU, typename NUM, typename DER, typename IVT, typename OVT> 407 friend int generic_gespmv_threaded (const SpMat<IU,NUM,DER> & A, const int32_t * indx, const IVT * numx, int32_t nnzx, 408 int32_t * & sendindbuf, OVT * & sendnumbuf, int * & sdispls, int p_c); 409 }; 410 411 // At this point, complete type of of SpDCCols is known, safe to declare these specialization (but macros won't work as they are preprocessed) 412 // General case #1: When both NT is the same 413 template <class IT, class NT> struct promote_trait< SpDCCols<IT,NT> , SpDCCols<IT,NT> > 414 { 415 typedef SpDCCols<IT,NT> T_promote; 416 }; 417 // General case #2: First is boolean the second is anything except boolean (to prevent ambiguity) 418 template <class IT, class NT> struct promote_trait< SpDCCols<IT,bool> , SpDCCols<IT,NT>, typename combblas::disable_if< combblas::is_boolean<NT>::value >::type > 419 { 420 typedef SpDCCols<IT,NT> T_promote; 421 }; 422 // General case #3: Second is boolean the first is anything except boolean (to prevent ambiguity) 423 template <class IT, class NT> struct promote_trait< SpDCCols<IT,NT> , SpDCCols<IT,bool>, typename combblas::disable_if< combblas::is_boolean<NT>::value >::type > 424 { 425 typedef SpDCCols<IT,NT> T_promote; 426 }; 427 template <class IT> struct promote_trait< SpDCCols<IT,int> , SpDCCols<IT,float> > 428 { 429 typedef SpDCCols<IT,float> T_promote; 430 }; 431 432 template <class IT> struct promote_trait< SpDCCols<IT,float> , SpDCCols<IT,int> > 433 { 434 typedef SpDCCols<IT,float> T_promote; 435 }; 436 template <class IT> struct promote_trait< SpDCCols<IT,int> , SpDCCols<IT,double> > 437 { 438 typedef SpDCCols<IT,double> T_promote; 439 }; 440 template <class IT> struct promote_trait< SpDCCols<IT,double> , SpDCCols<IT,int> > 441 { 442 typedef SpDCCols<IT,double> T_promote; 443 }; 444 445 446 // Below are necessary constructs to be able to define a SpMat<NT,IT> where 447 // all we know is DER (say SpDCCols<int, double>) and NT,IT 448 // in other words, we infer the templated SpDCCols<> type 449 // This is not a type conversion from an existing object, 450 // but a type inference for the newly created object 451 // NIT: New IT, NNT: New NT 452 template <class DER, class NIT, class NNT> 453 struct create_trait 454 { 455 // none 456 }; 457 458 // Capture everything of the form SpDCCols<OIT, ONT> 459 // it may come as a surprise that the partial specializations can 460 // involve more template parameters than the primary template 461 template <class NIT, class NNT, class OIT, class ONT> 462 struct create_trait< SpDCCols<OIT, ONT> , NIT, NNT > 463 { 464 typedef SpDCCols<NIT,NNT> T_inferred; 465 }; 466 467 } 468 469 #include "SpDCCols.cpp" 470 471 #endif 472