Epetra_MpiComm.cpp

Go to the documentation of this file.
00001 
00002 //@HEADER
00003 // ************************************************************************
00004 // 
00005 //               Epetra: Linear Algebra Services Package 
00006 //                 Copyright (2001) Sandia Corporation
00007 // 
00008 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00009 // license for use of this work by or on behalf of the U.S. Government.
00010 // 
00011 // This library is free software; you can redistribute it and/or modify
00012 // it under the terms of the GNU Lesser General Public License as
00013 // published by the Free Software Foundation; either version 2.1 of the
00014 // License, or (at your option) any later version.
00015 //  
00016 // This library is distributed in the hope that it will be useful, but
00017 // WITHOUT ANY WARRANTY; without even the implied warranty of
00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 // Lesser General Public License for more details.
00020 //  
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License along with this library; if not, write to the Free Software
00023 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 // USA
00025 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00026 // 
00027 // ************************************************************************
00028 //@HEADER
00029 
00030 #include "Epetra_MpiComm.h"
00031 
00032 //=============================================================================
00033 Epetra_MpiComm::Epetra_MpiComm(MPI_Comm Comm) : 
00034   Epetra_Object("Epetra::MpiComm"),
00035   MpiCommData_(new Epetra_MpiCommData(Comm))
00036 {
00037 }
00038 
00039 //=============================================================================
00040 Epetra_MpiComm::Epetra_MpiComm(const Epetra_MpiComm & Comm) : 
00041   Epetra_Object(Comm.Label()), 
00042   MpiCommData_(Comm.MpiCommData_)
00043 {
00044   MpiCommData_->IncrementReferenceCount();
00045 }
00046 
00047 //=============================================================================
00048 void Epetra_MpiComm::Barrier() const {
00049   MPI_Barrier(MpiCommData_->Comm_);
00050 }
00051 //=============================================================================
00052 int Epetra_MpiComm::Broadcast(double * Values, int Count, int Root) const {
00053   EPETRA_CHK_ERR(CheckInput(Values,Count));
00054   EPETRA_CHK_ERR(MPI_Bcast(Values, Count, MPI_DOUBLE, Root, MpiCommData_->Comm_));
00055   return(0);
00056 }
00057 //=============================================================================
00058 int Epetra_MpiComm::Broadcast(int * Values, int Count, int Root) const {
00059   EPETRA_CHK_ERR(CheckInput(Values,Count));
00060   EPETRA_CHK_ERR(MPI_Bcast(Values, Count, MPI_INT, Root, MpiCommData_->Comm_));
00061   return(0);
00062 }
00063 //=============================================================================
00064 int Epetra_MpiComm::Broadcast(long * Values, int Count, int Root) const {
00065   EPETRA_CHK_ERR(CheckInput(Values,Count));
00066   EPETRA_CHK_ERR(MPI_Bcast(Values, Count, MPI_LONG, Root, MpiCommData_->Comm_));
00067   return(0);
00068 }
00069 //=============================================================================
00070 int Epetra_MpiComm::GatherAll(double * MyVals, double * AllVals, int Count) const {
00071   EPETRA_CHK_ERR(CheckInput(MyVals,Count));
00072   EPETRA_CHK_ERR(CheckInput(AllVals,Count));
00073   EPETRA_CHK_ERR(MPI_Allgather(MyVals, Count, MPI_DOUBLE, AllVals, Count, MPI_DOUBLE, MpiCommData_->Comm_));
00074   return(0);
00075 }
00076 //=============================================================================
00077 int Epetra_MpiComm::GatherAll(int * MyVals, int * AllVals, int Count) const {
00078   EPETRA_CHK_ERR(CheckInput(MyVals,Count));
00079   EPETRA_CHK_ERR(CheckInput(AllVals,Count));
00080   EPETRA_CHK_ERR(MPI_Allgather(MyVals, Count, MPI_INT, AllVals, Count, MPI_INT, MpiCommData_->Comm_)); 
00081   return(0);
00082 }
00083 //=============================================================================
00084 int Epetra_MpiComm::GatherAll(long * MyVals, long * AllVals, int Count) const {
00085   EPETRA_CHK_ERR(CheckInput(MyVals,Count));
00086   EPETRA_CHK_ERR(CheckInput(AllVals,Count));
00087   EPETRA_CHK_ERR(MPI_Allgather(MyVals, Count, MPI_LONG, AllVals, Count, MPI_LONG, MpiCommData_->Comm_)); 
00088   return(0);
00089 }
00090 //=============================================================================
00091 int Epetra_MpiComm::SumAll(double * PartialSums, double * GlobalSums, int Count) const {
00092   EPETRA_CHK_ERR(CheckInput(PartialSums,Count));
00093   EPETRA_CHK_ERR(CheckInput(GlobalSums,Count));
00094   EPETRA_CHK_ERR(MPI_Allreduce(PartialSums, GlobalSums, Count, MPI_DOUBLE, MPI_SUM, MpiCommData_->Comm_));
00095   return(0);
00096 }
00097 //=============================================================================
00098 int Epetra_MpiComm::SumAll(int * PartialSums, int * GlobalSums, int Count) const {
00099   EPETRA_CHK_ERR(CheckInput(PartialSums,Count));
00100   EPETRA_CHK_ERR(CheckInput(GlobalSums,Count));
00101   EPETRA_CHK_ERR(MPI_Allreduce(PartialSums, GlobalSums, Count, MPI_INT, MPI_SUM, MpiCommData_->Comm_));
00102   return(0);
00103 }
00104 //=============================================================================
00105 int Epetra_MpiComm::SumAll(long * PartialSums, long * GlobalSums, int Count) const {
00106   EPETRA_CHK_ERR(CheckInput(PartialSums,Count));
00107   EPETRA_CHK_ERR(CheckInput(GlobalSums,Count));
00108   EPETRA_CHK_ERR(MPI_Allreduce(PartialSums, GlobalSums, Count, MPI_INT, MPI_SUM, MpiCommData_->Comm_));
00109   return(0);
00110 }
00111 //=============================================================================
00112 int Epetra_MpiComm::MaxAll(double * PartialMaxs, double * GlobalMaxs, int Count) const {
00113   EPETRA_CHK_ERR(CheckInput(PartialMaxs,Count));
00114   EPETRA_CHK_ERR(CheckInput(GlobalMaxs,Count));
00115   EPETRA_CHK_ERR(MPI_Allreduce(PartialMaxs, GlobalMaxs, Count, MPI_DOUBLE, MPI_MAX, MpiCommData_->Comm_));
00116   return(0);
00117 }
00118 //=============================================================================
00119 int Epetra_MpiComm::MaxAll(int * PartialMaxs, int * GlobalMaxs, int Count) const {
00120   EPETRA_CHK_ERR(CheckInput(PartialMaxs,Count));
00121   EPETRA_CHK_ERR(CheckInput(GlobalMaxs,Count));
00122   EPETRA_CHK_ERR(MPI_Allreduce(PartialMaxs, GlobalMaxs, Count, MPI_INT, MPI_MAX, MpiCommData_->Comm_));
00123   return(0);
00124 }
00125 //=============================================================================
00126 int Epetra_MpiComm::MaxAll(long * PartialMaxs, long * GlobalMaxs, int Count) const {
00127   EPETRA_CHK_ERR(CheckInput(PartialMaxs,Count));
00128   EPETRA_CHK_ERR(CheckInput(GlobalMaxs,Count));
00129   EPETRA_CHK_ERR(MPI_Allreduce(PartialMaxs, GlobalMaxs, Count, MPI_LONG, MPI_MAX, MpiCommData_->Comm_));
00130   return(0);
00131 }
00132 //=============================================================================
00133 int Epetra_MpiComm::MinAll(double * PartialMins, double * GlobalMins, int Count) const {
00134   EPETRA_CHK_ERR(CheckInput(PartialMins,Count));
00135   EPETRA_CHK_ERR(CheckInput(GlobalMins,Count));
00136   EPETRA_CHK_ERR(MPI_Allreduce(PartialMins, GlobalMins, Count, MPI_DOUBLE, MPI_MIN, MpiCommData_->Comm_));
00137   return(0);
00138 }
00139 //=============================================================================
00140 int Epetra_MpiComm::MinAll(int * PartialMins, int * GlobalMins, int Count) const {
00141   EPETRA_CHK_ERR(CheckInput(PartialMins,Count));
00142   EPETRA_CHK_ERR(CheckInput(GlobalMins,Count));
00143   EPETRA_CHK_ERR(MPI_Allreduce(PartialMins, GlobalMins, Count, MPI_INT, MPI_MIN, MpiCommData_->Comm_));
00144   return(0);
00145 }
00146 //=============================================================================
00147 int Epetra_MpiComm::MinAll(long * PartialMins, long * GlobalMins, int Count) const {
00148   EPETRA_CHK_ERR(CheckInput(PartialMins,Count));
00149   EPETRA_CHK_ERR(CheckInput(GlobalMins,Count));
00150   EPETRA_CHK_ERR(MPI_Allreduce(PartialMins, GlobalMins, Count, MPI_INT, MPI_MIN, MpiCommData_->Comm_));
00151   return(0);
00152 }
00153 //=============================================================================
00154 int Epetra_MpiComm::ScanSum(double * MyVals, double * ScanSums, int Count) const {
00155   EPETRA_CHK_ERR(CheckInput(MyVals,Count));
00156   EPETRA_CHK_ERR(CheckInput(ScanSums,Count));
00157   EPETRA_CHK_ERR(MPI_Scan(MyVals, ScanSums, Count, MPI_DOUBLE, MPI_SUM, MpiCommData_->Comm_));
00158   return(0);
00159 }
00160 //=============================================================================
00161 int Epetra_MpiComm::ScanSum(int * MyVals, int * ScanSums, int Count) const {
00162   EPETRA_CHK_ERR(CheckInput(MyVals,Count));
00163   EPETRA_CHK_ERR(CheckInput(ScanSums,Count));
00164   EPETRA_CHK_ERR(MPI_Scan(MyVals, ScanSums, Count, MPI_INT, MPI_SUM, MpiCommData_->Comm_));
00165   return(0);
00166 }
00167 //=============================================================================
00168 int Epetra_MpiComm::ScanSum(long * MyVals, long * ScanSums, int Count) const {
00169   EPETRA_CHK_ERR(CheckInput(MyVals,Count));
00170   EPETRA_CHK_ERR(CheckInput(ScanSums,Count));
00171   EPETRA_CHK_ERR(MPI_Scan(MyVals, ScanSums, Count, MPI_LONG, MPI_SUM, MpiCommData_->Comm_));
00172   return(0);
00173 }
00174 //=============================================================================
00175 Epetra_Distributor * Epetra_MpiComm:: CreateDistributor() const {
00176 
00177   Epetra_Distributor * dist = new Epetra_MpiDistributor(*this);
00178   return(dist);
00179 }
00180 //=============================================================================
00181 Epetra_Directory * Epetra_MpiComm:: CreateDirectory(const Epetra_BlockMap & map) const {
00182 
00183   Epetra_Directory * dir = new Epetra_BasicDirectory(map);
00184   return(dir);
00185 }
00186 //=============================================================================
00187 Epetra_MpiComm::~Epetra_MpiComm()  {
00188   CleanupData();
00189 }
00190 //=============================================================================
00191 void Epetra_MpiComm::CleanupData() {
00192   if(MpiCommData_ != 0) {
00193     MpiCommData_->DecrementReferenceCount();
00194     if(MpiCommData_->ReferenceCount() == 0) {
00195       delete MpiCommData_;
00196       MpiCommData_ = 0;
00197     }
00198   }
00199 }
00200 //=============================================================================
00201 Epetra_MpiComm & Epetra_MpiComm::operator= (const Epetra_MpiComm & Comm) {
00202   if((this != &Comm) && (MpiCommData_ != Comm.MpiCommData_)) {
00203     CleanupData();
00204     MpiCommData_ = Comm.MpiCommData_;
00205     MpiCommData_->IncrementReferenceCount();
00206   }
00207   return(*this);
00208 }

Generated on Thu Sep 18 12:37:58 2008 for Epetra Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1