00001 #ifndef ML_SPACE_H
00002 #define ML_SPACE_H
00003
00013
00014
00015
00016
00017
00018 #include "ml_common.h"
00019 #include "ml_comm.h"
00020 #include "MLAPI_Error.h"
00021 #include "MLAPI_Workspace.h"
00022 #include "MLAPI_BaseObject.h"
00023 #include "Teuchos_RefCountPtr.hpp"
00024 #include "Epetra_Map.h"
00025 #include "Epetra_IntSerialDenseVector.h"
00026 #include <iomanip>
00027
00028 namespace MLAPI {
00029
00040 class Space : public BaseObject {
00041
00042 public:
00043
00045
00047 Space()
00048 {
00049 NumMyElements_ = 0;
00050 NumGlobalElements_ = 0;
00051 IsLinear_ = false;
00052 Offset_ = 0;
00053 }
00054
00056
00066 Space(const int NumGlobalElements, const int NumMyElements = -1)
00067 {
00068 Reshape(NumGlobalElements, NumMyElements);
00069 }
00070
00072 Space(const Epetra_Map& Map)
00073 {
00074 Reshape(Map.NumGlobalElements(), Map.NumMyElements(),
00075 Map.MyGlobalElements());
00076 }
00077
00079
00088 Space(const int NumGlobalElements, const int NumMyElements,
00089 const int* MyGlobalElements)
00090 {
00091 Reshape(NumGlobalElements, NumMyElements, MyGlobalElements);
00092 }
00093
00095 Space(const Space& RHS)
00096 {
00097 NumMyElements_ = RHS.GetNumMyElements();
00098 NumGlobalElements_ = RHS.GetNumGlobalElements();
00099 Offset_ = RHS.GetOffset();
00100 IsLinear_ = RHS.IsLinear();
00101 RCPMyGlobalElements_ = RHS.GetRCPMyGlobalElements();
00102 }
00103
00105 ~Space() {};
00106
00108
00109
00111 void Reshape()
00112 {
00113 NumMyElements_ = 0;
00114 NumGlobalElements_ = 0;
00115 IsLinear_ = true;
00116 Offset_ = -1;
00117 RCPMyGlobalElements_ = Teuchos::null;
00118 }
00119
00121 void Reshape(const int NumGlobalElements, const int NumMyElements = -1)
00122 {
00123
00124 if (NumGlobalElements <= 0 && NumMyElements < 0)
00125 ML_THROW("NumGlobalElements = " + GetString(NumGlobalElements) +
00126 " and NumMyElements = " + GetString(NumMyElements), -1);
00127
00128 if (NumMyElements == -1) {
00129 NumMyElements_ = NumGlobalElements / GetNumProcs();
00130 if (GetMyPID() == 0)
00131 NumMyElements_ += NumGlobalElements % GetNumProcs();
00132 }
00133 else
00134 NumMyElements_ = NumMyElements;
00135
00136 NumGlobalElements_ = ML_Comm_GsumInt(GetML_Comm(),NumMyElements_);
00137
00138 if (NumGlobalElements != -1) {
00139 if (NumGlobalElements != NumGlobalElements_)
00140 ML_THROW("Specified # of global elements the sum of local elements (" +
00141 GetString(NumGlobalElements) + " vs. " +
00142 GetString(NumGlobalElements_), -1);
00143 }
00144
00145 Offset_ = ML_gpartialsum_int(NumMyElements_,GetML_Comm());
00146 IsLinear_ = true;
00147
00148 }
00149
00151 void Reshape(const int NumGlobalElements, const int NumMyElements,
00152 const int* MyGlobalElements)
00153 {
00154 if (NumGlobalElements <= 0 && NumMyElements < 0)
00155 ML_THROW("NumGlobalElements = " + GetString(NumGlobalElements) +
00156 " and NumMyElements = " + GetString(NumMyElements), -1);
00157
00158 if (NumMyElements == -1) {
00159 NumMyElements_ = NumGlobalElements / GetNumProcs();
00160 if (GetMyPID() == 0)
00161 NumMyElements_ += NumGlobalElements % GetNumProcs();
00162 }
00163 else
00164 NumMyElements_ = NumMyElements;
00165
00166 NumGlobalElements_ = ML_Comm_GsumInt(GetML_Comm(),NumMyElements_);
00167
00168 if (NumGlobalElements != -1) {
00169 if (NumGlobalElements != NumGlobalElements_)
00170 ML_THROW("Specified # of global elements the sum of local elements (" +
00171 GetString(NumGlobalElements) + " vs. " +
00172 GetString(NumGlobalElements_), -1);
00173 }
00174
00175 RCPMyGlobalElements_ = Teuchos::rcp(new Epetra_IntSerialDenseVector);
00176 RCPMyGlobalElements_->Resize(NumMyElements_);
00177 for (int i = 0 ; i < NumMyElements_ ; ++i)
00178 (*RCPMyGlobalElements_)[i] = MyGlobalElements[i];
00179
00180 Offset_ = -1;
00181 IsLinear_ = false;
00182 }
00183
00185
00186
00188 Space& operator=(const Space& RHS)
00189 {
00190 NumMyElements_ = RHS.GetNumMyElements();
00191 NumGlobalElements_ = RHS.GetNumGlobalElements();
00192 Offset_ = RHS.GetOffset();
00193 IsLinear_ = RHS.IsLinear();
00194 RCPMyGlobalElements_ = RHS.GetRCPMyGlobalElements();
00195 return(*this);
00196 }
00197
00199
00200
00201
00202 inline bool operator== (const Space& RHS) const
00203 {
00204 if (IsLinear() != RHS.IsLinear())
00205 return(false);
00206 if (GetNumGlobalElements() != RHS.GetNumGlobalElements())
00207 return(false);
00208 if (GetNumMyElements() != RHS.GetNumMyElements())
00209 return(false);
00210
00211 return(true);
00212 }
00213
00215
00216
00217
00218 inline bool operator!= (const Space& RHS) const
00219 {
00220 return(!(*this == RHS));
00221 }
00222
00224 Space& operator=(const string& Label)
00225 {
00226 SetLabel(Label);
00227 return(*this);
00228 }
00229
00231 inline int operator() (int i) const
00232 {
00233 #ifdef MLAPI_CHECK
00234 if ((i < 0) || (i >= NumMyElements_))
00235 ML_THROW("not valid index, " + GetString(i), -1);
00236 #endif
00237 if (IsLinear())
00238 return(i + GetOffset());
00239 else
00240 return((*RCPMyGlobalElements_)[i]);
00241 }
00242
00243
00244
00245
00247 int GetNumMyElements() const
00248 {
00249 return(NumMyElements_);
00250 }
00251
00253 int GetNumGlobalElements() const
00254 {
00255 return(NumGlobalElements_);
00256 }
00257
00259 inline int GetOffset() const
00260 {
00261 return(Offset_);
00262 }
00263
00265 inline bool IsLinear() const
00266 {
00267 return(IsLinear_);
00268 }
00269
00271 inline const Teuchos::RefCountPtr<Epetra_IntSerialDenseVector>
00272 GetRCPMyGlobalElements() const
00273 {
00274 return(RCPMyGlobalElements_);
00275 }
00276
00277
00278
00279
00281 std::ostream& Print(std::ostream& os,
00282 const bool verbose = true) const
00283 {
00284 if (GetMyPID() == 0) {
00285 os << endl;
00286 os << "*** MLAPI::Space ***" << endl;
00287 os << "Label = " << GetLabel() << endl;
00288 os << "NumMyElements() = " << GetNumMyElements() << endl;
00289 os << "NumGlobalElements() = " << GetNumGlobalElements() << endl;
00290 os << "Offset = " << GetOffset() << endl;
00291 if (IsLinear())
00292 os << "Distribution is linear" << endl;
00293 else
00294 os << "Distribution is not linear" << endl;
00295 os << endl;
00296 }
00297
00298 if (verbose) {
00299
00300 for (int iproc = 0 ; iproc < GetNumProcs() ; ++iproc) {
00301
00302 if (GetMyPID() == iproc) {
00303
00304 if (GetMyPID() == 0) {
00305 os.width(10);
00306 os << "ProcID";
00307 os.width(20);
00308 os << "LID";
00309 os.width(20);
00310 os << "GID" << endl << endl;
00311 }
00312
00313 for (int i = 0 ; i < GetNumMyElements() ; ++i) {
00314 os.width(10);
00315 os << GetMyPID();
00316 os.width(20);
00317 os << i;
00318 os.width(20);
00319 os << (*this)(i) << endl;
00320 }
00321 }
00322 Barrier();
00323 }
00324 Barrier();
00325
00326 if (GetMyPID() == 0)
00327 os << endl;
00328 }
00329
00330 Barrier();
00331
00332 return(os);
00333 }
00334
00335
00336
00337 private:
00338
00340 int NumMyElements_;
00342 int NumGlobalElements_;
00344 bool IsLinear_;
00346 int Offset_;
00348 Teuchos::RefCountPtr<Epetra_IntSerialDenseVector> RCPMyGlobalElements_;
00349
00350 };
00351
00352 }
00353
00354 #endif