44 #include "Epetra_Map.h" 51 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 54 const vector<int> & RowStencil,
58 BaseGraph_( BaseGraph ),
59 RowStencil_int_( vector< vector<int> >(1,RowStencil) ),
60 RowIndices_int_( vector<int>(1,rowIndex) ),
67 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 70 const vector<long long> & RowStencil,
74 BaseGraph_( BaseGraph ),
75 RowStencil_LL_( vector< vector<long long> >(1,RowStencil) ),
76 RowIndices_LL_( vector<long long>(1,rowIndex) ),
84 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 87 const vector< vector<int> > & RowStencil,
88 const vector<int> & RowIndices,
91 BaseGraph_( BaseGraph ),
92 RowStencil_int_( RowStencil ),
93 RowIndices_int_( RowIndices ),
99 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 102 const vector< vector<long long> > & RowStencil,
103 const vector<long long> & RowIndices,
106 BaseGraph_( BaseGraph ),
107 RowStencil_LL_( RowStencil ),
108 RowIndices_LL_( RowIndices ),
121 BaseGraph_( BaseGraph ),
122 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
126 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
130 ROffset_(
BlockUtility::CalculateOffset64(BaseGraph.RowMap())),
131 COffset_(
BlockUtility::CalculateOffset64(BaseGraph.ColMap()))
133 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 138 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 143 throw "EpetraExt::BlockCrsMatrix::BlockCrsMatrix: Error, Global indices unknown.";
147 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 150 const vector< vector<int> > & RowStencil,
151 const vector<int> & RowIndices,
154 BaseGraph_(
Copy, BaseMatrix.RowMatrixRowMap(), 1 ),
155 RowStencil_int_( RowStencil ),
156 RowIndices_int_( RowIndices ),
157 ROffset_(
BlockUtility::CalculateOffset64(BaseMatrix.RowMatrixRowMap())),
158 COffset_(
BlockUtility::CalculateOffset64(BaseMatrix.RowMatrixColMap()))
163 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 166 const vector< vector<long long> > & RowStencil,
167 const vector<long long> & RowIndices,
169 :
Epetra_CrsMatrix(
Copy, *(BlockUtility::GenerateBlockGraph( BaseMatrix, RowStencil, RowIndices, GlobalComm )) ),
170 BaseGraph_(
Copy, BaseMatrix.RowMatrixRowMap(), 1 ),
171 RowStencil_LL_( RowStencil ),
172 RowIndices_LL_( RowIndices ),
173 ROffset_(BlockUtility::CalculateOffset64(BaseMatrix.RowMatrixRowMap())),
174 COffset_(BlockUtility::CalculateOffset64(BaseMatrix.RowMatrixColMap()))
182 BaseGraph_( Matrix.BaseGraph_ ),
183 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
184 RowStencil_int_( Matrix.RowStencil_int_ ),
185 RowIndices_int_( Matrix.RowIndices_int_ ),
187 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
188 RowStencil_LL_( Matrix.RowStencil_LL_ ),
189 RowIndices_LL_( Matrix.RowIndices_LL_ ),
191 ROffset_( Matrix.ROffset_ ),
192 COffset_( Matrix.COffset_ )
202 template<
typename int_type>
203 void BlockCrsMatrix::TLoadBlock(
const Epetra_RowMatrix & BaseMatrix,
const int_type Row,
const int_type Col)
205 std::vector<int_type>& RowIndices_ = TRowIndices<int_type>();
206 std::vector< std::vector<int_type> >& RowStencil_ = TRowStencil<int_type>();
207 int_type RowOffset = RowIndices_[(std::size_t)Row] *
ROffset_;
208 int_type ColOffset = (RowIndices_[(std::size_t)Row] + RowStencil_[(std::size_t)Row][(std::size_t)Col]) *
COffset_;
219 vector<int> Indices_local(MaxIndices);
220 vector<int_type> Indices_global(MaxIndices);
221 vector<double> vals(MaxIndices);
226 BaseMatrix.
ExtractMyRowCopy( i, MaxIndices, NumIndices, &vals[0], &Indices_local[0] );
229 for(
int l = 0; l < NumIndices; ++l )
230 Indices_global[l] = ColOffset + (int_type) BaseColMap.GID64(Indices_local[l]);
232 int_type BaseRow = (int_type) BaseMap.GID64(i);
233 ierr = this->
ReplaceGlobalValues(BaseRow + RowOffset, NumIndices, &vals[0], &Indices_global[0]);
234 if (ierr != 0) std::cout <<
"WARNING BlockCrsMatrix::LoadBlock ReplaceGlobalValues err = " << ierr <<
235 "\n\t Row " << BaseRow + RowOffset <<
"Col start" << Indices_global[0] << std::endl;
240 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 244 return TLoadBlock<int>(BaseMatrix, Row, Col);
246 throw "EpetraExt::BlockCrsMatrix::LoadBlock: Global indices not int";
250 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 254 return TLoadBlock<long long>(BaseMatrix, Row, Col);
256 throw "EpetraExt::BlockCrsMatrix::LoadBlock: Global indices not long long";
261 template<
typename int_type>
262 void BlockCrsMatrix::TSumIntoBlock(
double alpha,
const Epetra_RowMatrix & BaseMatrix,
const int_type Row,
const int_type Col)
264 std::vector<int_type>& RowIndices_ = TRowIndices<int_type>();
265 std::vector< std::vector<int_type> >& RowStencil_ = TRowStencil<int_type>();
266 int_type RowOffset = RowIndices_[(std::size_t)Row] *
ROffset_;
267 int_type ColOffset = (RowIndices_[(std::size_t)Row] + RowStencil_[(std::size_t)Row][(std::size_t)Col]) *
COffset_;
278 vector<int> Indices_local(MaxIndices);
279 vector<int_type> Indices_global(MaxIndices);
280 vector<double> vals(MaxIndices);
285 BaseMatrix.
ExtractMyRowCopy( i, MaxIndices, NumIndices, &vals[0], &Indices_local[0] );
288 for(
int l = 0; l < NumIndices; ++l ) {
289 Indices_global[l] = ColOffset + (int_type) BaseColMap.GID64(Indices_local[l]);
293 int_type BaseRow = (int_type) BaseMap.GID64(i);
294 ierr = this->
SumIntoGlobalValues(BaseRow + RowOffset, NumIndices, &vals[0], &Indices_global[0]);
296 std::cout <<
"WARNING BlockCrsMatrix::SumIntoBlock SumIntoGlobalValues " 297 "err = " << ierr << std::endl <<
"\t Row " << BaseRow + RowOffset <<
298 "Col start" << Indices_global[0] << std::endl;
303 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 307 return TSumIntoBlock<int>(alpha, BaseMatrix, Row, Col);
309 throw "EpetraExt::BlockCrsMatrix::SumIntoBlock: Global indices not int";
313 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 317 return TSumIntoBlock<long long>(alpha, BaseMatrix, Row, Col);
319 throw "EpetraExt::BlockCrsMatrix::SumIntoBlock: Global indices not long long";
324 template<
typename int_type>
325 void BlockCrsMatrix::TSumIntoGlobalBlock(
double alpha,
const Epetra_RowMatrix & BaseMatrix,
const int_type Row,
const int_type Col)
327 int_type RowOffset = Row *
ROffset_;
328 int_type ColOffset = Col *
COffset_;
339 vector<int> Indices_local(MaxIndices);
340 vector<int_type> Indices_global(MaxIndices);
341 vector<double> vals(MaxIndices);
346 BaseMatrix.
ExtractMyRowCopy( i, MaxIndices, NumIndices, &vals[0], &Indices_local[0] );
349 for(
int l = 0; l < NumIndices; ++l ) {
350 Indices_global[l] = ColOffset + (int_type) BaseColMap.GID64(Indices_local[l]);
354 int_type BaseRow = (int_type) BaseMap.GID64(i);
355 ierr = this->
SumIntoGlobalValues(BaseRow + RowOffset, NumIndices, &vals[0], &Indices_global[0]);
357 std::cout <<
"WARNING BlockCrsMatrix::SumIntoBlock SumIntoGlobalValues " 358 <<
"err = " << ierr << std::endl
359 <<
"\t Row " << BaseRow + RowOffset
360 <<
" Col start" << Indices_global[0] << std::endl;
365 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 369 return TSumIntoGlobalBlock<int>(alpha, BaseMatrix, Row, Col);
371 throw "EpetraExt::BlockCrsMatrix::SumIntoGlobalBlock: Global indices not int";
375 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 379 return TSumIntoGlobalBlock<long long>(alpha, BaseMatrix, Row, Col);
381 throw "EpetraExt::BlockCrsMatrix::SumIntoGlobalBlock: Global indices not long long";
387 template<
typename int_type>
388 void BlockCrsMatrix::TBlockSumIntoGlobalValues(
const int_type BaseRow,
int NumIndices,
389 double* vals,
const int_type* Indices,
const int_type Row,
const int_type Col)
392 std::vector<int_type>& RowIndices_ = TRowIndices<int_type>();
393 std::vector< std::vector<int_type> >& RowStencil_ = TRowStencil<int_type>();
394 int_type RowOffset = RowIndices_[(std::size_t)Row] *
ROffset_;
395 int_type ColOffset = (RowIndices_[(std::size_t)Row] + RowStencil_[(std::size_t)Row][(std::size_t)Col]) *
COffset_;
398 vector<int_type> OffsetIndices(NumIndices);
399 for(
int l = 0; l < NumIndices; ++l ) OffsetIndices[l] = Indices[l] + ColOffset;
402 vals, &OffsetIndices[0]);
405 std::cout <<
"WARNING BlockCrsMatrix::BlockSumIntoGlobalValues err = " 406 << ierr << std::endl <<
"\t Row " << BaseRow + RowOffset
407 <<
" Col start" << Indices[0] << std::endl;
411 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 413 double* vals,
const int* Indices,
const int Row,
const int Col)
416 return TBlockSumIntoGlobalValues<int>(BaseRow, NumIndices, vals, Indices, Row, Col);
418 throw "EpetraExt::BlockCrsMatrix::BlockSumIntoGlobalValues: Global indices not int";
422 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 424 double* vals,
const long long* Indices,
const long long Row,
const long long Col)
427 return TBlockSumIntoGlobalValues<long long>(BaseRow, NumIndices, vals, Indices, Row, Col);
429 throw "EpetraExt::BlockCrsMatrix::BlockSumIntoGlobalValues: Global indices not long long";
434 template<
typename int_type>
435 void BlockCrsMatrix::TBlockReplaceGlobalValues(
const int_type BaseRow,
int NumIndices,
436 double* vals,
const int_type* Indices,
const int_type Row,
const int_type Col)
439 std::vector<int_type>& RowIndices_ = TRowIndices<int_type>();
440 std::vector< std::vector<int_type> >& RowStencil_ = TRowStencil<int_type>();
441 int_type RowOffset = RowIndices_[(std::size_t)Row] *
ROffset_;
442 int_type ColOffset = (RowIndices_[(std::size_t)Row] + RowStencil_[(std::size_t)Row][(std::size_t)Col]) *
COffset_;
445 vector<int_type> OffsetIndices(NumIndices);
446 for(
int l = 0; l < NumIndices; ++l ) OffsetIndices[l] = Indices[l] + ColOffset;
449 vals, &OffsetIndices[0]);
452 std::cout <<
"WARNING BlockCrsMatrix::BlockReplaceGlobalValues err = " 453 << ierr <<
"\n\t Row " << BaseRow + RowOffset <<
"Col start" 454 << Indices[0] << std::endl;
458 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 460 double* vals,
const int* Indices,
const int Row,
const int Col)
463 return TBlockReplaceGlobalValues<int>(BaseRow, NumIndices, vals, Indices, Row, Col);
465 throw "EpetraExt::BlockCrsMatrix::BlockReplaceGlobalValues: Global indices not int";
469 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 471 double* vals,
const long long* Indices,
const long long Row,
const long long Col)
474 return TBlockReplaceGlobalValues<long long>(BaseRow, NumIndices, vals, Indices, Row, Col);
476 throw "EpetraExt::BlockCrsMatrix::BlockReplaceGlobalValues: Global indices not long long";
481 template<
typename int_type>
482 void BlockCrsMatrix::TBlockExtractGlobalRowView(
const int_type BaseRow,
489 std::vector<int_type>& RowIndices_ = TRowIndices<int_type>();
490 std::vector< std::vector<int_type> >& RowStencil_ = TRowStencil<int_type>();
491 int_type RowOffset = RowIndices_[(std::size_t)Row] *
ROffset_;
492 int_type ColOffset = (RowIndices_[(std::size_t)Row] + RowStencil_[(std::size_t)Row][(std::size_t)Col]) *
COffset_;
500 NumEntries -= ColOffset;
503 std::cout <<
"WARNING BlockCrsMatrix::BlockExtractGlobalRowView err = " 504 << ierr <<
"\n\t Row " << BaseRow + RowOffset
505 <<
" Col " << Col+ColOffset << std::endl;
509 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 511 double*& vals,
const int Row,
const int Col)
514 return TBlockExtractGlobalRowView<int>(BaseRow, NumEntries, vals, Row, Col);
516 throw "EpetraExt::BlockCrsMatrix::BlockExtractGlobalRowView: Global indices not int";
520 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 522 double*& vals,
const long long Row,
const long long Col)
525 return TBlockExtractGlobalRowView<long long>(BaseRow, NumEntries, vals, Row, Col);
527 throw "EpetraExt::BlockCrsMatrix::BlockExtractGlobalRowView: Global indices not long long";
533 template<
typename int_type>
534 void BlockCrsMatrix::TExtractBlock(
Epetra_CrsMatrix & BaseMatrix,
const int_type Row,
const int_type Col)
536 std::vector<int_type>& RowIndices_ = TRowIndices<int_type>();
537 std::vector< std::vector<int_type> >& RowStencil_ = TRowStencil<int_type>();
538 int_type RowOffset = RowIndices_[(std::size_t)Row] *
ROffset_;
539 int_type ColOffset = (RowIndices_[(std::size_t)Row] + RowStencil_[(std::size_t)Row][(std::size_t)Col]) *
COffset_;
550 vector<int_type> Indices(MaxIndices);
551 vector<double> vals(MaxIndices);
563 int_type BaseRow = (int_type) BaseMap.GID64(i);
565 ierr = this->
ExtractMyRowView(myBlkBaseRow, BlkNumIndices, BlkValues, BlkIndices);
569 for(
int l = 0; l < BlkNumIndices; ++l ) {
571 indx = icol - ColOffset;
573 Indices[NumIndices] = indx;
574 vals[NumIndices] = BlkValues[l];
584 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 588 return TExtractBlock<int>(BaseMatrix, Row, Col);
590 throw "EpetraExt::BlockCrsMatrix::ExtractBlock: Global indices not int";
594 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 598 return TExtractBlock<long long>(BaseMatrix, Row, Col);
600 throw "EpetraExt::BlockCrsMatrix::ExtractBlock: Global indices not long long";
int MaxNumEntries() const
void LoadBlock(const Epetra_RowMatrix &BaseMatrix, const int Row, const int Col)
Routine for loading a base matrices values into the large Block Matrix The Row and Col arguments are ...
virtual int ReplaceGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
virtual const Epetra_Map & RowMatrixRowMap() const=0
virtual int SumIntoGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
int ExtractGlobalRowView(int GlobalRow, int &NumEntries, double *&Values, int *&Indices) const
std::vector< std::vector< int > > RowStencil_int_
static void GenerateRowStencil(const Epetra_CrsGraph &LocalBlockGraph, std::vector< int > RowIndices, std::vector< std::vector< int > > &RowStencil)
Generate stencil arrays from a local block graph.
void BlockReplaceGlobalValues(const int BaseRow, int NumIndices, double *Values, const int *Indices, const int Row, const int Col)
std::vector< long long > RowIndices_LL_
void ExtractBlock(Epetra_CrsMatrix &BaseMatrix, const int Row, const int Col)
const Epetra_BlockMap & RowMap() const
int ExtractMyRowView(int MyRow, int &NumEntries, double *&Values, int *&Indices) const
EpetraExt::BlockCrsMatrix: A class for constructing a distributed block matrix.
const Epetra_Map & RowMatrixRowMap() const
static long long CalculateOffset64(const Epetra_BlockMap &BaseMap)
void BlockExtractGlobalRowView(const int BaseRow, int &NumEntries, double *&Values, const int Row, const int Col)
const Epetra_BlockMap & ColMap() const
virtual int MaxNumEntries() const=0
std::vector< std::vector< long long > > RowStencil_LL_
std::vector< int > RowIndices_int_
bool GlobalIndicesInt() const
void SumIntoGlobalBlock(double alpha, const Epetra_RowMatrix &BaseMatrix, const int Row, const int Col)
Routine for summing base matrices values into the large Block Matrix The Row and Col arguments are gl...
int NumMyElements() const
virtual ~BlockCrsMatrix()
Destructor.
const Epetra_Map & RowMatrixColMap() const
BlockCrsMatrix(const Epetra_CrsGraph &BaseGraph, const std::vector< int > &RowStencil, int RowIndex, const Epetra_Comm &GlobalComm)
BlockCrsMatrix constuctor with one block row per processor.
void SumIntoBlock(double alpha, const Epetra_RowMatrix &BaseMatrix, const int Row, const int Col)
Routine for summing base matrices values into the large Block Matrix The Row and Col arguments are in...
virtual int ExtractMyRowCopy(int MyRow, int Length, int &NumEntries, double *Values, int *Indices) const=0
void BlockSumIntoGlobalValues(const int BaseRow, int NumIndices, double *Values, const int *Indices, const int Row, const int Col)
Sum Entries into Block matrix using base-matrix numbering plus block Row and Col The Row and Col argu...
virtual const Epetra_Map & RowMatrixColMap() const=0
bool GlobalIndicesLongLong() const
static Epetra_CrsGraph * GenerateBlockGraph(const Epetra_CrsGraph &BaseGraph, const std::vector< std::vector< int > > &RowStencil, const std::vector< int > &RowIndices, const Epetra_Comm &GlobalComm)
BlockCrsMatrix constuctor.