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_ )
202template<
typename int_type>
203void 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_;
218 int MaxIndices = BaseMatrix.MaxNumEntries();
219 vector<int> Indices_local(MaxIndices);
220 vector<int_type> Indices_global(MaxIndices);
221 vector<double> vals(MaxIndices);
225 for (
int i=0; i<BaseMap.NumMyElements(); i++) {
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";
261template<
typename int_type>
262void 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_;
277 int MaxIndices = BaseMatrix.MaxNumEntries();
278 vector<int> Indices_local(MaxIndices);
279 vector<int_type> Indices_global(MaxIndices);
280 vector<double> vals(MaxIndices);
284 for (
int i=0; i<BaseMap.NumMyElements(); i++) {
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";
324template<
typename int_type>
325void 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_;
338 int MaxIndices = BaseMatrix.MaxNumEntries();
339 vector<int> Indices_local(MaxIndices);
340 vector<int_type> Indices_global(MaxIndices);
341 vector<double> vals(MaxIndices);
345 for (
int i=0; i<BaseMap.NumMyElements(); i++) {
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";
387template<
typename int_type>
388void 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";
434template<
typename int_type>
435void 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";
481template<
typename int_type>
482void 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";
533template<
typename int_type>
534void 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_;
549 int MaxIndices = BaseMatrix.MaxNumEntries();
550 vector<int_type> Indices(MaxIndices);
551 vector<double> vals(MaxIndices);
560 for (
int i=0; i<BaseMap.NumMyElements(); i++) {
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];
580 BaseMatrix.ReplaceGlobalValues(BaseRow, NumIndices, &vals[0], &Indices[0]);
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";