EpetraExt Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
EpetraExt_RestrictedMultiVectorWrapper.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 "EpetraExt_ConfigDefs.h"
45 
46 
47 #ifdef HAVE_MPI
48 
49 
51 #include "Epetra_MpiComm.h"
52 #include "Epetra_BlockMap.h"
53 #include "Epetra_MultiVector.h"
54 
55 
56 namespace EpetraExt {
57 
58 
59 RestrictedMultiVectorWrapper::RestrictedMultiVectorWrapper()
60  : proc_is_active(true),
61  subcomm_is_set(false),
62  MPI_SubComm_(MPI_COMM_NULL),
63  RestrictedComm_(0),
64  ResMap_(0),
65  input_mv_(),
66  restricted_mv_()
67  {}
68 
69 
70 RestrictedMultiVectorWrapper::~RestrictedMultiVectorWrapper() {
71  delete ResMap_;
72  delete RestrictedComm_;
73 }
74 
75 
76 int RestrictedMultiVectorWrapper::SetMPISubComm(MPI_Comm MPI_SubComm){
77  if(!subcomm_is_set){
78  MPI_SubComm_=MPI_SubComm; delete RestrictedComm_; subcomm_is_set=true;
79  return 0;
80  }
81  else return -1;
82 }
83 
84 
85 int
86 RestrictedMultiVectorWrapper::
87 restrict_comm (Teuchos::RCP<Epetra_MultiVector> input_mv)
88 {
89  using Teuchos::rcp;
90  input_mv_ = input_mv;
91 
92  // Extract the input MV's communicator and Map.
93  const Epetra_MpiComm *InComm =
94  dynamic_cast<const Epetra_MpiComm*> (& (input_mv_->Comm ()));
95  const Epetra_BlockMap *InMap =
96  dynamic_cast<const Epetra_BlockMap*> (& (input_mv_->Map ()));
97 
98  if (! InComm || ! InMap) {
99  return -1; // At least one dynamic cast failed.
100  }
101 
102  if (! subcomm_is_set) {
103  /* Build the Split Communicators, If Needed */
104  int color;
105  if (InMap->NumMyElements()) {
106  color = 1;
107  }
108  else {
109  color = MPI_UNDEFINED;
110  }
111  const int err =
112  MPI_Comm_split (InComm->Comm(), color, InComm->MyPID(), &MPI_SubComm_);
113  if (err != MPI_SUCCESS) {
114  return -2;
115  }
116  }
117  else {
118  /* Sanity check user-provided subcomm - drop an error if the MPISubComm
119  does not include a processor with data. */
120  if (input_mv->MyLength() && MPI_SubComm_ == MPI_COMM_NULL) {
121  return -2;
122  }
123  }
124 
125  /* Mark active processors */
126  if (MPI_SubComm_ == MPI_COMM_NULL) {
127  proc_is_active = false;
128  }
129  else {
130  proc_is_active = true;
131  }
132 
133  if (proc_is_active) {
134 
135 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
136  if(InMap->GlobalIndicesInt()) {
137  int Nrows = InMap->NumGlobalElements ();
138  RestrictedComm_ = new Epetra_MpiComm (MPI_SubComm_);
139 
140  // Build the restricted Maps
141  ResMap_ = new Epetra_BlockMap (Nrows, InMap->NumMyElements(),
142  InMap->MyGlobalElements(),
143  InMap->ElementSizeList(),
144  InMap->IndexBase(), *RestrictedComm_);
145  }
146  else
147 #endif
148 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
149  if(InMap->GlobalIndicesLongLong()) {
150  long long Nrows = InMap->NumGlobalElements64 ();
151  RestrictedComm_ = new Epetra_MpiComm (MPI_SubComm_);
152 
153  // Build the restricted Maps
154  ResMap_ = new Epetra_BlockMap (Nrows, InMap->NumMyElements(),
155  InMap->MyGlobalElements64(),
156  InMap->ElementSizeList(),
157  InMap->IndexBase64(), *RestrictedComm_);
158  }
159  else
160 #endif
161  throw "EpetraExt::RestrictedMultiVectorWrapper::restrict_comm ERROR: GlobalIndices type unknown";
162 
163  // Allocate the restricted matrix
164  double *A;
165  int LDA;
166  input_mv_->ExtractView (&A,&LDA);
167  restricted_mv_ = rcp (new Epetra_MultiVector (View, *ResMap_, A, LDA,
168  input_mv_->NumVectors ()));
169  }
170  return 0; // Success!
171 }/*end restrict_comm*/
172 
173 }
174 
175 #endif
int NumGlobalElements() const
bool GlobalIndicesLongLong() const
long long IndexBase64() const
int MyGlobalElements(int *MyGlobalElementList) const
long long NumGlobalElements64() const
int * ElementSizeList() const
bool GlobalIndicesInt() const
int IndexBase() const
int NumMyElements() const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
MPI_Comm Comm() const
long long * MyGlobalElements64() const
int MyPID() const