Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Epetra_MpiComm.cpp
Go to the documentation of this file.
1 
2 //@HEADER
3 // ************************************************************************
4 //
5 // Epetra: 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 #include "Epetra_MpiComm.h"
44 
45 //=============================================================================
46 Epetra_MpiComm::Epetra_MpiComm(MPI_Comm theComm) :
47  Epetra_Object("Epetra::MpiComm"),
48  MpiCommData_(new Epetra_MpiCommData(theComm))
49 {
50 }
51 
52 //=============================================================================
54  Epetra_Object(theComm.Label()),
55  MpiCommData_(theComm.MpiCommData_)
56 {
58 }
59 
60 //=============================================================================
62  MPI_Barrier(MpiCommData_->Comm_);
63 }
64 //=============================================================================
65 int Epetra_MpiComm::Broadcast(double * Values, int Count, int Root) const {
66  EPETRA_CHK_ERR(CheckInput(Values,Count));
67  EPETRA_CHK_ERR(MPI_Bcast(Values, Count, MPI_DOUBLE, Root, MpiCommData_->Comm_));
68  return(0);
69 }
70 //=============================================================================
71 int Epetra_MpiComm::Broadcast(int * Values, int Count, int Root) const {
72  EPETRA_CHK_ERR(CheckInput(Values,Count));
73  EPETRA_CHK_ERR(MPI_Bcast(Values, Count, MPI_INT, Root, MpiCommData_->Comm_));
74  return(0);
75 }
76 //=============================================================================
77 int Epetra_MpiComm::Broadcast(long * Values, int Count, int Root) const {
78  EPETRA_CHK_ERR(CheckInput(Values,Count));
79  EPETRA_CHK_ERR(MPI_Bcast(Values, Count, MPI_LONG, Root, MpiCommData_->Comm_));
80  return(0);
81 }
82 //=============================================================================
83 int Epetra_MpiComm::Broadcast(long long * Values, int Count, int Root) const {
84  EPETRA_CHK_ERR(CheckInput(Values,Count));
85  EPETRA_CHK_ERR(MPI_Bcast(Values, Count, MPI_LONG_LONG, Root, MpiCommData_->Comm_));
86  return(0);
87 }
88 //=============================================================================
89 int Epetra_MpiComm::Broadcast(char * Values, int Count, int Root) const {
90  EPETRA_CHK_ERR(CheckInput(Values,Count));
91  EPETRA_CHK_ERR(MPI_Bcast(Values, Count, MPI_CHAR, Root, MpiCommData_->Comm_));
92  return(0);
93 }
94 //=============================================================================
95 int Epetra_MpiComm::GatherAll(double * MyVals, double * AllVals, int Count) const {
96  EPETRA_CHK_ERR(CheckInput(MyVals,Count));
97  EPETRA_CHK_ERR(CheckInput(AllVals,Count));
98  EPETRA_CHK_ERR(MPI_Allgather(MyVals, Count, MPI_DOUBLE, AllVals, Count, MPI_DOUBLE, MpiCommData_->Comm_));
99  return(0);
100 }
101 //=============================================================================
102 int Epetra_MpiComm::GatherAll(int * MyVals, int * AllVals, int Count) const {
103  EPETRA_CHK_ERR(CheckInput(MyVals,Count));
104  EPETRA_CHK_ERR(CheckInput(AllVals,Count));
105  EPETRA_CHK_ERR(MPI_Allgather(MyVals, Count, MPI_INT, AllVals, Count, MPI_INT, MpiCommData_->Comm_));
106  return(0);
107 }
108 //=============================================================================
109 int Epetra_MpiComm::GatherAll(long * MyVals, long * AllVals, int Count) const {
110  EPETRA_CHK_ERR(CheckInput(MyVals,Count));
111  EPETRA_CHK_ERR(CheckInput(AllVals,Count));
112  EPETRA_CHK_ERR(MPI_Allgather(MyVals, Count, MPI_LONG, AllVals, Count, MPI_LONG, MpiCommData_->Comm_));
113  return(0);
114 }
115 //=============================================================================
116 int Epetra_MpiComm::GatherAll(long long * MyVals, long long * AllVals, int Count) const {
117  EPETRA_CHK_ERR(CheckInput(MyVals,Count));
118  EPETRA_CHK_ERR(CheckInput(AllVals,Count));
119  EPETRA_CHK_ERR(MPI_Allgather(MyVals, Count, MPI_LONG_LONG, AllVals, Count, MPI_LONG_LONG, MpiCommData_->Comm_));
120  return(0);
121 }
122 //=============================================================================
123 int Epetra_MpiComm::SumAll(double * PartialSums, double * GlobalSums, int Count) const {
124  EPETRA_CHK_ERR(CheckInput(PartialSums,Count));
125  EPETRA_CHK_ERR(CheckInput(GlobalSums,Count));
126  EPETRA_CHK_ERR(MPI_Allreduce(PartialSums, GlobalSums, Count, MPI_DOUBLE, MPI_SUM, MpiCommData_->Comm_));
127  return(0);
128 }
129 //=============================================================================
130 int Epetra_MpiComm::SumAll(int * PartialSums, int * GlobalSums, int Count) const {
131  EPETRA_CHK_ERR(CheckInput(PartialSums,Count));
132  EPETRA_CHK_ERR(CheckInput(GlobalSums,Count));
133  EPETRA_CHK_ERR(MPI_Allreduce(PartialSums, GlobalSums, Count, MPI_INT, MPI_SUM, MpiCommData_->Comm_));
134  return(0);
135 }
136 //=============================================================================
137 int Epetra_MpiComm::SumAll(long * PartialSums, long * GlobalSums, int Count) const {
138  EPETRA_CHK_ERR(CheckInput(PartialSums,Count));
139  EPETRA_CHK_ERR(CheckInput(GlobalSums,Count));
140  EPETRA_CHK_ERR(MPI_Allreduce(PartialSums, GlobalSums, Count, MPI_LONG, MPI_SUM, MpiCommData_->Comm_));
141  return(0);
142 }
143 //=============================================================================
144 int Epetra_MpiComm::SumAll(long long * PartialSums, long long * GlobalSums, int Count) const {
145  EPETRA_CHK_ERR(CheckInput(PartialSums,Count));
146  EPETRA_CHK_ERR(CheckInput(GlobalSums,Count));
147  EPETRA_CHK_ERR(MPI_Allreduce(PartialSums, GlobalSums, Count, MPI_LONG_LONG, MPI_SUM, MpiCommData_->Comm_));
148  return(0);
149 }
150 //=============================================================================
151 int Epetra_MpiComm::MaxAll(double * PartialMaxs, double * GlobalMaxs, int Count) const {
152  EPETRA_CHK_ERR(CheckInput(PartialMaxs,Count));
153  EPETRA_CHK_ERR(CheckInput(GlobalMaxs,Count));
154  EPETRA_CHK_ERR(MPI_Allreduce(PartialMaxs, GlobalMaxs, Count, MPI_DOUBLE, MPI_MAX, MpiCommData_->Comm_));
155  return(0);
156 }
157 //=============================================================================
158 int Epetra_MpiComm::MaxAll(int * PartialMaxs, int * GlobalMaxs, int Count) const {
159  EPETRA_CHK_ERR(CheckInput(PartialMaxs,Count));
160  EPETRA_CHK_ERR(CheckInput(GlobalMaxs,Count));
161  EPETRA_CHK_ERR(MPI_Allreduce(PartialMaxs, GlobalMaxs, Count, MPI_INT, MPI_MAX, MpiCommData_->Comm_));
162  return(0);
163 }
164 //=============================================================================
165 int Epetra_MpiComm::MaxAll(long * PartialMaxs, long * GlobalMaxs, int Count) const {
166  EPETRA_CHK_ERR(CheckInput(PartialMaxs,Count));
167  EPETRA_CHK_ERR(CheckInput(GlobalMaxs,Count));
168  EPETRA_CHK_ERR(MPI_Allreduce(PartialMaxs, GlobalMaxs, Count, MPI_LONG, MPI_MAX, MpiCommData_->Comm_));
169  return(0);
170 }
171 //=============================================================================
172 int Epetra_MpiComm::MaxAll(long long * PartialMaxs, long long * GlobalMaxs, int Count) const {
173  EPETRA_CHK_ERR(CheckInput(PartialMaxs,Count));
174  EPETRA_CHK_ERR(CheckInput(GlobalMaxs,Count));
175  EPETRA_CHK_ERR(MPI_Allreduce(PartialMaxs, GlobalMaxs, Count, MPI_LONG_LONG, MPI_MAX, MpiCommData_->Comm_));
176  return(0);
177 }
178 //=============================================================================
179 int Epetra_MpiComm::MinAll(double * PartialMins, double * GlobalMins, int Count) const {
180  EPETRA_CHK_ERR(CheckInput(PartialMins,Count));
181  EPETRA_CHK_ERR(CheckInput(GlobalMins,Count));
182  EPETRA_CHK_ERR(MPI_Allreduce(PartialMins, GlobalMins, Count, MPI_DOUBLE, MPI_MIN, MpiCommData_->Comm_));
183  return(0);
184 }
185 //=============================================================================
186 int Epetra_MpiComm::MinAll(int * PartialMins, int * GlobalMins, int Count) const {
187  EPETRA_CHK_ERR(CheckInput(PartialMins,Count));
188  EPETRA_CHK_ERR(CheckInput(GlobalMins,Count));
189  EPETRA_CHK_ERR(MPI_Allreduce(PartialMins, GlobalMins, Count, MPI_INT, MPI_MIN, MpiCommData_->Comm_));
190  return(0);
191 }
192 //=============================================================================
193 int Epetra_MpiComm::MinAll(long * PartialMins, long * GlobalMins, int Count) const {
194  EPETRA_CHK_ERR(CheckInput(PartialMins,Count));
195  EPETRA_CHK_ERR(CheckInput(GlobalMins,Count));
196  EPETRA_CHK_ERR(MPI_Allreduce(PartialMins, GlobalMins, Count, MPI_LONG, MPI_MIN, MpiCommData_->Comm_));
197  return(0);
198 }
199 //=============================================================================
200 int Epetra_MpiComm::MinAll(long long * PartialMins, long long * GlobalMins, int Count) const {
201  EPETRA_CHK_ERR(CheckInput(PartialMins,Count));
202  EPETRA_CHK_ERR(CheckInput(GlobalMins,Count));
203  EPETRA_CHK_ERR(MPI_Allreduce(PartialMins, GlobalMins, Count, MPI_LONG_LONG, MPI_MIN, MpiCommData_->Comm_));
204  return(0);
205 }
206 //=============================================================================
207 int Epetra_MpiComm::ScanSum(double * MyVals, double * ScanSums, int Count) const {
208  EPETRA_CHK_ERR(CheckInput(MyVals,Count));
209  EPETRA_CHK_ERR(CheckInput(ScanSums,Count));
210  EPETRA_CHK_ERR(MPI_Scan(MyVals, ScanSums, Count, MPI_DOUBLE, MPI_SUM, MpiCommData_->Comm_));
211  return(0);
212 }
213 //=============================================================================
214 int Epetra_MpiComm::ScanSum(int * MyVals, int * ScanSums, int Count) const {
215  EPETRA_CHK_ERR(CheckInput(MyVals,Count));
216  EPETRA_CHK_ERR(CheckInput(ScanSums,Count));
217  EPETRA_CHK_ERR(MPI_Scan(MyVals, ScanSums, Count, MPI_INT, MPI_SUM, MpiCommData_->Comm_));
218  return(0);
219 }
220 //=============================================================================
221 int Epetra_MpiComm::ScanSum(long * MyVals, long * ScanSums, int Count) const {
222  EPETRA_CHK_ERR(CheckInput(MyVals,Count));
223  EPETRA_CHK_ERR(CheckInput(ScanSums,Count));
224  EPETRA_CHK_ERR(MPI_Scan(MyVals, ScanSums, Count, MPI_LONG, MPI_SUM, MpiCommData_->Comm_));
225  return(0);
226 }
227 //=============================================================================
228 int Epetra_MpiComm::ScanSum(long long * MyVals, long long * ScanSums, int Count) const {
229  EPETRA_CHK_ERR(CheckInput(MyVals,Count));
230  EPETRA_CHK_ERR(CheckInput(ScanSums,Count));
231  EPETRA_CHK_ERR(MPI_Scan(MyVals, ScanSums, Count, MPI_LONG_LONG, MPI_SUM, MpiCommData_->Comm_));
232  return(0);
233 }
234 //=============================================================================
236 
237  Epetra_Distributor * dist = new Epetra_MpiDistributor(*this);
238  return(dist);
239 }
240 //=============================================================================
242 
243  Epetra_Directory * dir = new Epetra_BasicDirectory(map);
244  return(dir);
245 }
246 //=============================================================================
248  CleanupData();
249 }
250 //=============================================================================
252  if(MpiCommData_ != 0) {
254  if(MpiCommData_->ReferenceCount() == 0) {
255  delete MpiCommData_;
256  MpiCommData_ = 0;
257  }
258  }
259 }
260 //=============================================================================
262  if((this != &theComm) && (MpiCommData_ != theComm.MpiCommData_)) {
263  CleanupData();
264  MpiCommData_ = theComm.MpiCommData_;
266  }
267  return(*this);
268 }
int CheckInput(double *ptr, int count) const
int GatherAll(double *MyVals, double *AllVals, int Count) const
Epetra_MpiComm All Gather function.
Epetra_Distributor: The Epetra Gather/Scatter Setup Base Class.
void DecrementReferenceCount()
Decrement reference count.
Definition: Epetra_Data.cpp:66
void Barrier() const
Epetra_MpiComm Barrier function.
#define EPETRA_CHK_ERR(a)
int MaxAll(double *PartialMaxs, double *GlobalMaxs, int Count) const
Epetra_MpiComm Global Max function.
int MinAll(double *PartialMins, double *GlobalMins, int Count) const
Epetra_MpiComm Global Min function.
Epetra_MpiComm & operator=(const Epetra_MpiComm &Comm)
Assignment Operator.
MPI implementation of Epetra_Distributor.
Epetra_MpiComm: The Epetra MPI Communication Class.
int ScanSum(double *MyVals, double *ScanSums, int Count) const
Epetra_MpiComm Scan Sum function.
Epetra_Distributor * CreateDistributor() const
Create a distributor object.
Epetra_Directory: This class is a pure virtual class whose interface allows Epetra_Map and Epetr_Bloc...
int Broadcast(double *MyVals, int Count, int Root) const
Epetra_MpiComm Broadcast function.
Epetra_BasicDirectory: This class allows Epetra_Map objects to reference non-local elements...
Epetra_Object: The base Epetra class.
Definition: Epetra_Object.h:57
int SumAll(double *PartialSums, double *GlobalSums, int Count) const
Epetra_MpiComm Global Sum function.
Epetra_MpiCommData * MpiCommData_
virtual ~Epetra_MpiComm()
Epetra_MpiComm Destructor.
int ReferenceCount() const
Get reference count.
Definition: Epetra_Data.cpp:71
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
void IncrementReferenceCount()
Increment reference count.
Definition: Epetra_Data.cpp:61
Epetra_Directory * CreateDirectory(const Epetra_BlockMap &Map) const
Create a directory object for the given Epetra_BlockMap.
Epetra_MpiCommData: The Epetra Mpi Communication Data Class.
Epetra_MpiComm(MPI_Comm comm)
Epetra_MpiComm MPI Constructor.