00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029 #include "EpetraExt_BlockMultiVector.h"
00030 #include "EpetraExt_BlockUtility.h"
00031 #include "Epetra_Map.h"
00032
00033 using std::vector;
00034
00035 namespace EpetraExt {
00036
00037
00038
00039 BlockMultiVector::BlockMultiVector(
00040 const Epetra_BlockMap & BaseMap,
00041 const Epetra_BlockMap & GlobalMap,
00042 int NumVectors )
00043 : Epetra_MultiVector( GlobalMap, NumVectors ),
00044 BaseMap_( BaseMap ),
00045 Offset_( BlockUtility::CalculateOffset( BaseMap ) )
00046 {
00047 }
00048
00049
00050
00051 BlockMultiVector::BlockMultiVector(
00052 Epetra_DataAccess CV,
00053 const Epetra_BlockMap & BaseMap,
00054 const Epetra_MultiVector & BlockVec)
00055 : Epetra_MultiVector( CV, BlockVec, 0, BlockVec.NumVectors() ),
00056 BaseMap_( BaseMap ),
00057 Offset_( BlockUtility::CalculateOffset( BaseMap ) )
00058 {
00059 }
00060
00061
00062
00063 BlockMultiVector::BlockMultiVector(const BlockMultiVector& Source)
00064 : Epetra_MultiVector( dynamic_cast<const Epetra_MultiVector &>(Source) ),
00065 BaseMap_( Source.BaseMap_ ),
00066 Offset_( Source.Offset_ )
00067 {
00068 }
00069
00070
00071 BlockMultiVector::~BlockMultiVector()
00072 {
00073 }
00074
00075
00076 int BlockMultiVector::ExtractBlockValues(Epetra_MultiVector & BaseVector, int GlobalBlockRow) const
00077 {
00078 int IndexOffset = GlobalBlockRow * Offset_;
00079 int localIndex=0;
00080
00081
00082
00083 for (int i=0; i<BaseMap_.NumMyElements(); i++) {
00084 localIndex = this->Map().LID((IndexOffset + BaseMap_.GID(i)));
00085 if (localIndex==-1) {
00086 cout << "Error in BlockMultiVector::GetBlock: " << i << " "
00087 << IndexOffset << " " << BaseMap_.GID(i) << endl;
00088 return -1;
00089 }
00090 for (int j=0; j<NumVectors(); j++)
00091 BaseVector[j][i] = (*this)[j][localIndex];
00092 }
00093
00094 return 0;
00095 }
00096
00097
00098 int BlockMultiVector::LoadBlockValues(const Epetra_MultiVector & BaseVector, int GlobalBlockRow)
00099 {
00100 int IndexOffset = GlobalBlockRow * Offset_;
00101 int localIndex=0;
00102
00103
00104
00105 for (int i=0; i<BaseMap_.NumMyElements(); i++) {
00106 localIndex = this->Map().LID((IndexOffset + BaseMap_.GID(i)));
00107 if (localIndex==-1) {
00108 cout << "Error in BlockMultiVector::GetBlock: " << i << " "
00109 << IndexOffset << " " << BaseMap_.GID(i) << endl;
00110 return -1;
00111 }
00112 for (int j=0; j<NumVectors(); j++)
00113 (*this)[j][localIndex] = BaseVector[j][i];
00114 }
00115
00116 return 0;
00117 }
00118
00119
00120 Teuchos::RCP<const Epetra_MultiVector>
00121 BlockMultiVector::GetBlock(int GlobalBlockRow) const
00122 {
00123 int offset = GlobalBlockRow * BaseMap_.NumMyElements();
00124 int numVecs = NumVectors();
00125 double **pointers = Pointers();
00126 double **block_pointers = new double*[numVecs];
00127 for (int i=0; i<numVecs; i++)
00128 block_pointers[i] = pointers[i]+offset;
00129 Teuchos::RCP<Epetra_MultiVector> block =
00130 Teuchos::rcp(new Epetra_MultiVector(View, BaseMap_, block_pointers,
00131 numVecs));
00132 delete [] block_pointers;
00133 return block;
00134 }
00135
00136
00137 Teuchos::RCP<Epetra_MultiVector>
00138 BlockMultiVector::GetBlock(int GlobalBlockRow)
00139 {
00140 int offset = GlobalBlockRow * BaseMap_.NumMyElements();
00141 int numVecs = NumVectors();
00142 double **pointers = Pointers();
00143 double **block_pointers = new double*[numVecs];
00144 for (int i=0; i<numVecs; i++)
00145 block_pointers[i] = pointers[i]+offset;
00146 Teuchos::RCP<Epetra_MultiVector> block =
00147 Teuchos::rcp(new Epetra_MultiVector(View, BaseMap_, block_pointers,
00148 numVecs));
00149 delete [] block_pointers;
00150 return block;
00151 }
00152
00153
00154 const Epetra_BlockMap&
00155 BlockMultiVector::GetBaseMap() const
00156 {
00157 return BaseMap_;
00158 }
00159
00160 }