EpetraExt Package Browser (Single Doxygen Collection)  Development
EpetraExt_PointToBlockDiagPermute.cpp
Go to the documentation of this file.
1 /*
2 //@HEADER
3 // ***********************************************************************
4 //
5 // EpetraExt: Epetra Extended - Linear Algebra Services Package
6 // Copyright (2011) Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ***********************************************************************
41 //@HEADER
42 */
43 
44 #include "Epetra_ConfigDefs.h"
47 #include "Epetra_MultiVector.h"
48 #include "Epetra_Import.h"
49 #include "Epetra_Export.h"
50 #include "Epetra_Comm.h"
51 
52 #include <stdio.h>
53 #include <fstream>
54 
55 #define MAX(x,y) ((x)>(y)?(x):(y))
56 
57 //=========================================================================
59  :Epetra_DistObject(MAT.RowMap()),
60  Matrix_(&MAT),
61  PurelyLocalMode_(true),
62  ContiguousBlockMode_(false),
63  ContiguousBlockSize_(0),
64  NumBlocks_(0),
65  Blockstart_(0),
66  Blockids_int_(0),
67 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
68  Blockids_LL_(0),
69 #endif
70  BDMap_(0),
71  CompatibleMap_(0),
72  BDMat_(0),
73  Importer_(0),
74  Exporter_(0),
75  ImportVector_(0),
76  ExportVector_(0)
77 {
78 
79 }
80 
81 //=========================================================================
82 // Destructor
84 {
85  if(BDMap_) delete BDMap_;
86  if(CompatibleMap_) delete CompatibleMap_;
87  if(BDMat_) delete BDMat_;
88  if(Importer_) delete Importer_;
89  if(Exporter_) delete Exporter_;
90  if(ImportVector_) delete ImportVector_;
91  if(ExportVector_) delete ExportVector_;
92 }
93 
94 //=========================================================================
95 // Set list
96 template<typename int_type>
98  List_=List;
99 
100  // Check for contiguous blocking first
101  ContiguousBlockSize_=List_.get("contiguous block size",0);
102  if(ContiguousBlockSize_!=0){
104  PurelyLocalMode_=false;
105  }
106 
107  // Local vs. global ids & mode
108  NumBlocks_=List_.get("number of local blocks",0);
109  Blockstart_=List_.get("block start index",(int*)0);
110  Blockids_ref<int>()=List_.get("block entry lids",(int*)0);
111  if(Blockids_const_ptr<int>())
112  PurelyLocalMode_=true;
113  else{
114  Blockids_ref<int_type>()=List_.get("block entry gids",(int_type*)0);
115  if(Blockids_const_ptr<int_type>()) {
116  PurelyLocalMode_=false;
117  if(sizeof(int) != sizeof(int_type)) {
118  Blockids_ref<int>() = new int[Matrix_->RowMap().NumMyElements()];
119  }
120  }
121  }
122 
123  // Sanity checks
125  // Can't use contiguous at the same time as the other modes
126  if(NumBlocks_ || Blockstart_ || Blockids_const_ptr<int>() || Blockids_const_ptr<int_type>()) EPETRA_CHK_ERR(-4);
127  }
128  else {
129  if(NumBlocks_ <= 0) EPETRA_CHK_ERR(-1);
130  if(!Blockstart_) EPETRA_CHK_ERR(-2);
131  if(!Blockids_const_ptr<int>() && !Blockids_const_ptr<int_type>()) EPETRA_CHK_ERR(-3);
132  }
133 
134  return 0;
135 }
136 
138 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
139  if(Matrix_->RowMap().GlobalIndicesInt()) {
140  return TSetParameters<int>(List);
141  }
142  else
143 #endif
144 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
146  return TSetParameters<long long>(List);
147  }
148  else
149 #endif
150  throw "EpetraExt_PointToBlockDiagPermute::SetParameters: ERROR, GlobalIndices type unknown.";
151 }
152 
153 //=========================================================================
154 // Extracts the block-diagonal, builds maps, etc.
156  int rv=ExtractBlockDiagonal();
157  return rv;
158 }
159 
160 //=========================================================================
161 // Returns the result of a Epetra_Operator applied to a Epetra_MultiVector X in Y.
163  return -1;
164 
165 }
166 
167 //=========================================================================
168 // Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector X in Y.
170  // Stuff borrowed from Epetra_CrsMatrix
171  int NumVectors = X.NumVectors();
172  if (NumVectors!=Y.NumVectors()) {
173  EPETRA_CHK_ERR(-2); // Need same number of vectors in each MV
174  }
175 
176  const Epetra_MultiVector *Xp=&X;
177  Epetra_MultiVector *Yp=&Y;
178 
179  // Allocate temp workspace if X==Y and there are no imports or exports
180  Epetra_MultiVector * Xcopy = 0;
181  if (&X==&Y && Importer_==0 && Exporter_==0) {
182  Xcopy = new Epetra_MultiVector(X);
183  Xp=Xcopy;
184  }
185 
186  UpdateImportVector(NumVectors); // Make sure Import and Export Vectors are compatible
187  UpdateExportVector(NumVectors);
188 
189  // If we have a non-trivial importer, we must import elements that are permuted or are on other processors
190  if (Importer_){
192  Xp=ImportVector_;
193  }
194 
195  // If we have a non-trivial exporter, we must export elements that are permuted or belong to other processors
196  if (Exporter_) {
197  Yp=ExportVector_;
198  }
199 
200  // Do the matvec
201  BDMat_->ApplyInverse(*Xp,*Yp);
202 
203  // Export if needed
204  if (Exporter_) {
205  Y.PutScalar(0.0); // Make sure target is zero
206  Y.Export(*ExportVector_, *Exporter_, Add); // Fill Y with Values from export vector
207  }
208 
209  // Cleanup
210  if(Xcopy) {
211  delete Xcopy;
212  EPETRA_CHK_ERR(1); // Return positive code to alert the user about needing extra copy of X
213  return 1;
214  }
215 
216  return 0;
217 }
218 
219 //=========================================================================
220 // Print method
221 void EpetraExt_PointToBlockDiagPermute::Print(std::ostream& os) const{
222  if(Importer_) std::cout<<*Importer_<<std::endl;
223  if(Exporter_) std::cout<<*Exporter_<<std::endl;
224  if(BDMat_) std::cout<<*BDMat_<<std::endl;
225 }
226 
227 //=========================================================================
228 // Pulls the block diagonal of the matrix and then builds the BDMat_
229 template<typename int_type>
231  int i,j;
232  std::vector<int_type> ExtRows;
233  int* LocalColIDS=0;
234  int ExtSize=0,ExtCols=0,MainCols=0;
235  int Nrows=Matrix_->NumMyRows();
236  int *l2b,*block_offset;
237  int_type *l2blockid;
238  const Epetra_Map &RowMap=Matrix_->RowMap();
239  int index,col,row_in_block,col_in_block,length,*colind;
240  double *values;
241 
242  bool verbose=(bool)(List_.get("output",0) > 0);
243 
244  // Contiguous Setup
246 
247  // Compute block size lists
248  int *bsize=new int[NumBlocks_];
249  for(i=0;i<NumBlocks_;i++)
250  bsize[i]=Blockstart_[i+1]-Blockstart_[i];
251 
252  // Use the ScanSum function to compute a prefix sum of the number of points
253  int_type MyMinGID;
254  int_type NumBlocks_int_type = NumBlocks_;
255  Matrix_->Comm().ScanSum(&NumBlocks_int_type,&MyMinGID, 1);
256  MyMinGID-=NumBlocks_;
257  int_type *MyBlockGIDs=new int_type[NumBlocks_];
258  for(i=0;i<NumBlocks_;i++)
259  MyBlockGIDs[i]=MyMinGID+i;
260 
261  BDMap_=new Epetra_BlockMap((int_type) -1,NumBlocks_,MyBlockGIDs,bsize,0,Matrix_->Comm());
262 
263  // Allocate the Epetra_DistBlockMatrix
265 
266  // First check to see if we can switch back to PurelyLocalMode_ if we're not
267  // in it
268  if(!PurelyLocalMode_){
269  // Find the non-local row IDs
270  ExtSize=MAX(0,Blockstart_[NumBlocks_]-RowMap.NumMyElements());// estimate
271  ExtRows.reserve(ExtSize);
272 
273  for(i=0;i<Blockstart_[NumBlocks_];i++){
274  if(RowMap.LID(Blockids_const_ptr<int_type>()[i])==-1)
275  ExtRows.push_back(Blockids_const_ptr<int_type>()[i]);
276  }
277 
278  // Switch to PurelyLocalMode_ if we never need GIDs
279  int size_sum;
280  ExtSize=ExtRows.size();
281  RowMap.Comm().SumAll(&ExtSize,&size_sum,1);
282  if(size_sum==0){
283  if(verbose && !Matrix_->Comm().MyPID()) printf("EpetraExt_PointToBlockDiagPermute: Switching to purely local mode\n");
284  PurelyLocalMode_=true;
285  for(i=0;i<Blockstart_[NumBlocks_];i++){
286  Blockids_ref<int>()[i]=RowMap.LID(Blockids_const_ptr<int_type>()[i]);
287  }
288  }
289  }
290 
291  if(PurelyLocalMode_){
292  /*******************************************************/
293  // Allocations
294  l2b=new int[Nrows];
295  block_offset=new int[Nrows];
296 
297  // Build the local-id-to-block-id list, block offset list
298  for(i=0;i<NumBlocks_;i++) {
299  for(j=Blockstart_[i];j<Blockstart_[i+1];j++){
300  block_offset[Blockids_const_ptr<int>()[j]]=j-Blockstart_[i];
301  l2b[Blockids_const_ptr<int>()[j]]=i;
302  }
303  }
304 
305  // Copy the data to the EpetraExt_BlockDiagMatrix
306  for(i=0;i<Nrows;i++) {
307  int block_num = l2b[i];
308  if(block_num>=0 && block_num<NumBlocks_) {
309  row_in_block=block_offset[i];
310  Matrix_->ExtractMyRowView(i,length,values,colind);
311  int Nnz=0;
312  for (j = 0; j < length; j++) {
313  col = colind[j];
314  if(col < Nrows && l2b[col]==block_num){
315  Nnz++;
316  col_in_block = block_offset[col];
317  index = col_in_block * bsize[block_num] + row_in_block;
318  (*BDMat_)[block_num][index]=values[j];
319  }
320  }
321  /* Handle the case of a zero row. */
322  /* By just putting a 1 on the diagonal. */
323  if (Nnz == 0) {
324  index = row_in_block * bsize[block_num] + row_in_block;
325  (*BDMat_)[block_num][index]=1.0;
326  }
327  }
328  }
329 
330  // Build the compatible map for import/export
331  l2blockid=new int_type[Nrows];
332  for(i=0;i<Nrows;i++)
333  l2blockid[Blockstart_[l2b[i]]+block_offset[i]]= (int_type) Matrix_->RowMap().GID64(i);
334 
335  // Build the Compatible Map, Import/Export Objects
336  CompatibleMap_=new Epetra_Map(-1,Nrows,l2blockid,0,Matrix_->Comm());
337 
338  }
339  else{
340  /*******************************************************/
341  // Do the import to grab matrix entries
342  // Allocate temporaries for import
343  int_type* ExtRowsPtr = ExtRows.size()>0 ? &ExtRows[0] : NULL;
344  Epetra_Map TmpMap((int_type) -1,ExtSize, ExtRowsPtr,0,Matrix_->Comm());;
345  Epetra_CrsMatrix TmpMatrix(Copy,TmpMap,0);
346  Epetra_Import TmpImporter(TmpMap,RowMap);
347 
348  TmpMatrix.Import(*Matrix_,TmpImporter,Insert);
349  TmpMatrix.FillComplete();
350 
351  ExtCols=TmpMatrix.NumMyCols();
352  MainCols=Matrix_->NumMyCols();
353 
354 
355  // Build the column reidex - main matrix
356  LocalColIDS=new int[MainCols+ExtCols];
357  for(i=0;i<MainCols;i++){
358  int_type GID= (int_type) Matrix_->GCID64(i);
359  int MainLID=RowMap.LID(GID);
360  if(MainLID!=-1) LocalColIDS[i]=MainLID;
361  else{
362  int ExtLID=TmpMatrix.LRID(GID);
363  if(ExtLID!=-1) LocalColIDS[i]=Nrows+ExtLID;
364  else LocalColIDS[i]=-1;
365  }
366  }
367 
368  // Build the column reidex - ext matrix
369  for(i=0;i<ExtCols;i++){
370  int_type GID= (int_type) TmpMatrix.GCID64(i);
371  int MainLID=RowMap.LID(GID);
372  if(MainLID!=-1) LocalColIDS[MainCols+i]=MainLID;
373  else{
374  int ExtLID=TmpMatrix.LRID(GID);
375  if(ExtLID!=-1) LocalColIDS[MainCols+i]=Nrows+ExtLID;
376  else LocalColIDS[MainCols+i]=-1;
377  }
378  }
379 
380  // Allocations
381  l2b=new int[Nrows+ExtSize];
382  block_offset=new int[Nrows+ExtSize];
383 
384  // Build l2b/block_offset with the expanded local index
385  //NTS: Uses the same ordering of operation as the above routine, which is why
386  //it works.
387  for(i=0;i<Nrows+ExtSize;i++) block_offset[i]=l2b[i]=-1;
388  int ext_idx=0;
389  for(i=0;i<NumBlocks_;i++) {
390  for(j=Blockstart_[i];j<Blockstart_[i+1];j++){
391  int LID=RowMap.LID(Blockids_const_ptr<int_type>()[j]);
392  if(LID==-1) {LID=Nrows+ext_idx;ext_idx++;}
393  block_offset[LID]=j-Blockstart_[i];
394  l2b[LID]=i;
395  }
396  }
397 
398  // Copy the data to the EpetraExt_BlockDiagMatrix from Matrix_
399  for(i=0;i<Nrows;i++) {
400  int block_num = l2b[i];
401  if(block_num>=0 && block_num<NumBlocks_) {
402  row_in_block=block_offset[i];
403  Matrix_->ExtractMyRowView(i,length,values,colind);
404  int Nnz=0;
405  for (j = 0; j < length; j++) {
406  col = LocalColIDS[colind[j]];
407  if(col!=-1 && l2b[col]==block_num){
408  Nnz++;
409  col_in_block = block_offset[col];
410  index = col_in_block * bsize[block_num] + row_in_block;
411  (*BDMat_)[block_num][index]=values[j];
412  }
413  }
414  /* Handle the case of a zero row. */
415  /* By just putting a 1 on the diagonal. */
416  if (Nnz == 0) {
417  index = row_in_block * bsize[block_num] + row_in_block;
418  (*BDMat_)[block_num][index]=values[j];
419  }
420  }
421  }
422 
423  // Copy the data to the EpetraExt_BlockDiagMatrix from TmpMatrix
424  for(i=0;i<ExtSize;i++) {
425  int block_num = l2b[Nrows+i];
426  if(block_num>=0 && block_num<NumBlocks_) {
427  row_in_block=block_offset[Nrows+i];
428  TmpMatrix.ExtractMyRowView(i,length,values,colind);
429  int Nnz=0;
430  for (j = 0; j < length; j++) {
431  col = LocalColIDS[MainCols+colind[j]];
432  if(col!=-1 && l2b[col]==block_num){
433  Nnz++;
434  col_in_block = block_offset[col];
435  index = col_in_block * bsize[block_num] + row_in_block;
436  (*BDMat_)[block_num][index]=values[j];
437  }
438  }
439  /* Handle the case of a zero row. */
440  /* By just putting a 1 on the diagonal. */
441  if (Nnz == 0) {
442  index = row_in_block * bsize[block_num] + row_in_block;
443  (*BDMat_)[block_num][index]=1.0;
444  }
445  }
446  }
447 
448  // Build the compatible map for import/export
449  l2blockid=new int_type[Blockstart_[NumBlocks_]];
450  for(i=0;i<Nrows;i++){
451  int bkid=l2b[i];
452  if(bkid>-1) l2blockid[Blockstart_[bkid]+block_offset[i]]= (int_type) RowMap.GID64(i);
453  }
454  // NTS: This is easier - if we imported it, we know we need it
455  for(i=0;i<ExtSize;i++)
456  l2blockid[Blockstart_[l2b[Nrows+i]]+block_offset[Nrows+i]]= (int_type) TmpMatrix.GRID64(i);
457 
458  // Build the Compatible Map, Import/Export Objects
460 
461  }//end else
462 
463  // Set BDMat_'s parameter list and compute!
464  Teuchos::ParameterList dummy,inList;
465  inList=List_.get("blockdiagmatrix: list",dummy);
466  BDMat_->SetParameters(inList);
467  BDMat_->Compute();
468 
469  // Build importer/exporter
474 
475  // Cleanup
476  delete [] LocalColIDS;
477  delete [] block_offset;
478  delete [] l2b;
479  delete [] l2blockid;
480  delete [] bsize;
481  delete [] MyBlockGIDs;
482 
483  // Contiguous Cleanup
485 
486  return 0;
487 }
488 
490 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
491  if(Matrix_->RowMap().GlobalIndicesInt()) {
492  return TExtractBlockDiagonal<int>();
493  }
494  else
495 #endif
496 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
498  return TExtractBlockDiagonal<long long>();
499  }
500  else
501 #endif
502  throw "EpetraExt_PointToBlockDiagPermute::ExtractBlockDiagonal: ERROR, GlobalIndices type unknown.";
503 }
504 
505 //=======================================================================================================
506 template<typename int_type>
508  if(!ContiguousBlockMode_) return 0;
509  // NTS: In case of processor-crossing blocks, the lowest PID always gets the block;
510  const Epetra_Map &RowMap=Matrix_->RowMap();
511 
512  int_type MinMyGID= (int_type) RowMap.MinMyGID64();
513  int_type MaxMyGID= (int_type) RowMap.MaxMyGID64();
514  int_type Base= (int_type) Matrix_->IndexBase64();
515 
516  // Find the GID that begins my first block
517  int_type MyFirstBlockGID=ContiguousBlockSize_*(int_type)ceil(((double)(MinMyGID - Base))/ContiguousBlockSize_)+Base;
518  NumBlocks_=(int)ceil((double)((MaxMyGID-MyFirstBlockGID+1.0)) / ContiguousBlockSize_);
519 
520  // Allocate memory
521  Blockstart_=new int[NumBlocks_+1];
522  Blockids_ref<int_type>()=new int_type[NumBlocks_*ContiguousBlockSize_];
523  if(sizeof(int) != sizeof(int_type))
524  Blockids_ref<int>()=new int[NumBlocks_*ContiguousBlockSize_];
526 
527  // Fill the arrays
528  for(int i=0,ct=0;i<NumBlocks_;i++){
529  Blockstart_[i]=ct;
530  for(int j=0;j<ContiguousBlockSize_;j++,ct++){
531  Blockids_ref<int_type>()[ct]=MyFirstBlockGID+ct;
532  }
533  }
534 
535  return 0;
536 }
537 
539 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
540  if(Matrix_->RowMap().GlobalIndicesInt()) {
541  return TSetupContiguousMode<int>();
542  }
543  else
544 #endif
545 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
547  return TSetupContiguousMode<long long>();
548  }
549  else
550 #endif
551  throw "EpetraExt_PointToBlockDiagPermute::SetupContiguousMode: ERROR, GlobalIndices type unknown.";
552 }
553 
554 //=======================================================================================================
556  if(!ContiguousBlockMode_) return 0;
557  NumBlocks_=0;
558  if(Blockstart_) {delete [] Blockstart_; Blockstart_=0;}
559  if(Blockids_int_) {delete [] Blockids_int_; Blockids_int_=0;}
560 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
561  if(Blockids_LL_) {delete [] Blockids_LL_; Blockids_LL_=0;}
562 #endif
563  return 0;
564 }
565 
566 
567 //=======================================================================================================
568 // Creates an Epetra_CrsMatrix from the BlockDiagMatrix. This is generally only useful if you want to do a matrix-matrix multiply.
569 template<typename int_type>
572 
573  const Epetra_BlockMap &BlockMap=BDMat_->BlockMap();
574  const Epetra_BlockMap &DataMap=BDMat_->DataMap();
575  const int *vlist=DataMap.FirstPointInElementList();
576  const int *xlist=BlockMap.FirstPointInElementList();
577  const int *blocksize=BlockMap.ElementSizeList();
578  const double *values=BDMat_->Values();
579  int NumBlocks=BDMat_->NumMyBlocks();
580 
581  // Maximum size vector for stashing GIDs
582  std::vector<int_type> GIDs;
583  GIDs.resize(BlockMap.MaxMyElementSize());
584 
585 
586  for(int i=0;i<NumBlocks;i++){
587  int Nb=blocksize[i];
588  int vidx0=vlist[i];
589  int xidx0=xlist[i];
590  // Get global indices
591  for(int j=0;j<Nb;j++)
592  GIDs[j]= (int_type) CompatibleMap_->GID64(xidx0+j);
593 
594  // Remember: We're using column-major storage for LAPACK's benefit
595  int ierr=NewMat->InsertGlobalValues(Nb,&GIDs[0],&values[vidx0],Epetra_FECrsMatrix::COLUMN_MAJOR);
596  if(ierr < 0) throw "CreateFECrsMatrix: ERROR in InsertGlobalValues";
597  }
598  NewMat->GlobalAssemble();
599  return NewMat;
600 }
601 
603 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
604  if(Matrix_->RowMap().GlobalIndicesInt()) {
605  return TCreateFECrsMatrix<int>();
606  }
607  else
608 #endif
609 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
611  return TCreateFECrsMatrix<long long>();
612  }
613  else
614 #endif
615  throw "EpetraExt_PointToBlockDiagPermute::CreateFECrsMatrix: ERROR, GlobalIndices type unknown.";
616 }
617 
618 //=======================================================================================================
620  if(Importer_ != 0) {
621  if(ImportVector_ != 0) {
622  if(ImportVector_->NumVectors() != NumVectors) {
623  delete ImportVector_;
624  ImportVector_= 0;
625  }
626  }
627  if(ImportVector_ == 0)
628  ImportVector_ = new Epetra_MultiVector(*CompatibleMap_,NumVectors); // Create import vector if needed
629  }
630  return;
631 }
632 //=======================================================================================================
634  if(Exporter_ != 0) {
635  if(ExportVector_ != 0) {
636  if(ExportVector_->NumVectors() != NumVectors) {
637  delete ExportVector_;
638  ExportVector_= 0;
639  }
640  }
641  if(ExportVector_ == 0)
642  ExportVector_ = new Epetra_MultiVector(*CompatibleMap_,NumVectors); // Create Export vector if needed
643  }
644  return;
645 }
646 
647 
648 
649 
650 //=========================================================================
652  return -1;
653 }
654 
655 //=========================================================================
657  return -1;
658 }
659 
660 //=========================================================================
662  return -1;
663 }
664 
665 //=========================================================================
667  return -1;
668 }
669 
670 //=========================================================================
671 // Allows the source and target (\e this) objects to be compared for compatibility, return nonzero if not.
673  return -1;
674 }
675 
676 //=========================================================================
677 // Perform ID copies and permutations that are on processor.
679  int NumSameIDs,
680  int NumPermuteIDs,
681  int * PermuteToLIDs,
682  int * PermuteFromLIDs,
683  const Epetra_OffsetIndex * Indexor,
684  Epetra_CombineMode CombineMode){
685  return -1;
686 }
687 
688 //=========================================================================
689 // Perform any packing or preparation required for call to DoTransfer().
691  int NumExportIDs,
692  int* ExportLIDs,
693  int& LenExports,
694  char*& Exports,
695  int& SizeOfPacket,
696  int* Sizes,
697  bool & VarSizes,
698  Epetra_Distributor& Distor){
699  return -1;
700 }
701 
702 //=========================================================================
703 // Perform any unpacking and combining after call to DoTransfer().
705  int NumImportIDs,
706  int* ImportLIDs,
707  int LenImports,
708  char* Imports,
709  int& SizeOfPacket,
710  Epetra_Distributor& Distor,
711  Epetra_CombineMode CombineMode,
712  const Epetra_OffsetIndex * Indexor){
713  return -1;
714 }
715 
716 
const Epetra_Map & RowMap() const
virtual void Print(std::ostream &os) const
Print information about this object to the given output stream.
int TSetParameters(Teuchos::ParameterList &List)
int NumMyBlocks() const
Returns the number of local blocks.
virtual int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector X in Y...
virtual int SetParameters(Teuchos::ParameterList &List)
Sets the parameter list.
long long MinMyGID64() const
#define MAX(x, y)
int LRID(int GRID_in) const
int NumMyRows() const
bool SameAs(const Epetra_BlockMap &Map) const
T & get(ParameterList &l, const std::string &name)
int NumMyCols() const
int PutScalar(double ScalarConstant)
#define EPETRA_CHK_ERR(a)
virtual const Epetra_BlockMap & DataMap() const
Returns the Epetra_BlockMap object with the distribution of underlying values.
Insert
int * FirstPointInElementList() const
const Epetra_Comm & Comm() const
int ExtractMyRowView(int MyRow, int &NumEntries, double *&Values, int *&Indices) const
virtual int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Epetra_Operator applied to a Epetra_MultiVector X in Y.
virtual int MyPID() const=0
const Epetra_Map & RangeMap() const
int FillComplete(bool OptimizeDataStorage=true)
int InsertGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
EpetraExt_BlockDiagMatrix: A class for storing distributed block matrices.
virtual int UnpackAndCombine(const Epetra_SrcDistObject &Source, int NumImportIDs, int *ImportLIDs, int LenImports, char *Imports, int &SizeOfPacket, Epetra_Distributor &Distor, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor)
Perform any unpacking and combining after call to DoTransfer().
virtual int SumAll(double *PartialSums, double *GlobalSums, int Count) const=0
int Export(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
long long GID64(int LID) const
virtual int CheckSizes(const Epetra_SrcDistObject &Source)
Allows the source and target (this) objects to be compared for compatibility, return nonzero if not...
int NumVectors() const
bool GlobalIndicesInt() const
int NumMyElements() const
virtual int Compute()
Computes the inverse / factorization if such is set on the list.
int GlobalAssemble(bool callFillComplete=true, Epetra_CombineMode combineMode=Add, bool save_off_and_reuse_map_exporter=false)
virtual int CopyAndPermute(const Epetra_SrcDistObject &Source, int NumSameIDs, int NumPermuteIDs, int *PermuteToLIDs, int *PermuteFromLIDs, const Epetra_OffsetIndex *Indexor, Epetra_CombineMode CombineMode=Zero)
Perform ID copies and permutations that are on processor.
const Epetra_Comm & Comm() const
long long GRID64(int LRID_in) const
virtual int SetParameters(Teuchos::ParameterList &List)
SetParameters.
int LID(int GID) const
EpetraExt_PointToBlockDiagPermute(const Epetra_CrsMatrix &MAT)
@ Name Constructors
virtual const Epetra_BlockMap & BlockMap() const
Returns the Epetra_BlockMap object associated with the range of this operator.
virtual int Compute()
Extracts the block-diagonal, builds maps, etc.
virtual int PackAndPrepare(const Epetra_SrcDistObject &Source, int NumExportIDs, int *ExportLIDs, int &LenExports, char *&Exports, int &SizeOfPacket, int *Sizes, bool &VarSizes, Epetra_Distributor &Distor)
Perform any packing or preparation required for call to DoTransfer().
virtual Epetra_FECrsMatrix * CreateFECrsMatrix()
Create an Epetra_FECrsMatrix from the BlockDiagMatrix.
int SetupContiguousMode()
Setup/Cleanup for Contiguous Mode.
Copy
Epetra_CombineMode
double * Values() const
Returns a pointer to the array containing the blocks.
long long MaxMyGID64() const
virtual int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector X in Y...
long long GCID64(int LCID_in) const
bool GlobalIndicesLongLong() const
int ExtractBlockDiagonal()
Pulls the block diagonal of the matrix and then builds the BDMat_.
Add
int Import(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
long long IndexBase64() const
virtual int ScanSum(double *MyVals, double *ScanSums, int Count) const=0
int Import(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Imports an Epetra_DistObject using the Epetra_Import object.
int Export(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Exports an Epetra_DistObject using the Epetra_Import object.
const Epetra_Map & DomainMap() const