FEI Version of the Day
fei_Graph_Impl.cpp
00001 /*
00002 // @HEADER
00003 // ************************************************************************
00004 //             FEI: Finite Element Interface to Linear Solvers
00005 //                  Copyright (2005) Sandia Corporation.
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the
00008 // U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Alan Williams (william@sandia.gov) 
00038 //
00039 // ************************************************************************
00040 // @HEADER
00041 */
00042 
00043 
00044 #include <fei_macros.hpp>
00045 
00046 #include <fei_Graph_Impl.hpp>
00047 #include <fei_EqnComm.hpp>
00048 #include <fei_CommUtils.hpp>
00049 #include <fei_TemplateUtils.hpp>
00050 #include <fei_VectorSpace.hpp>
00051 
00052 #undef fei_file
00053 #define fei_file "fei_Graph_Impl.cpp"
00054 #include <fei_ErrMacros.hpp>
00055 
00056 //----------------------------------------------------------------------------
00057 fei::Graph_Impl::Graph_Impl(MPI_Comm comm, int firstLocalRow, int lastLocalRow)
00058   : localGraphData_(NULL),
00059     remoteGraphData_(),
00060     eqnComm_(),
00061     firstLocalRow_(firstLocalRow),
00062     lastLocalRow_(lastLocalRow),
00063     localProc_(0),
00064     numProcs_(1),
00065     comm_(comm)
00066 {
00067   localProc_ = fei::localProc(comm_);
00068   numProcs_  = fei::numProcs(comm_);
00069   //for remoteGraphData_, we don't know what the range of row-numbers will
00070   //be, so we'll just construct it with -1,-1
00071   remoteGraphData_.resize(numProcs_);
00072   for(int p=0; p<numProcs_; ++p) {
00073     remoteGraphData_[p] = new remote_table_type(-1, -1);
00074   }
00075   eqnComm_.reset(new fei::EqnComm(comm_, lastLocalRow-firstLocalRow+1));
00076   localGraphData_       = new table_type(firstLocalRow_, lastLocalRow_);
00077 }
00078 
00079 //----------------------------------------------------------------------------
00080 fei::Graph_Impl::~Graph_Impl()
00081 {
00082   delete localGraphData_;
00083   for(unsigned i=0; i<remoteGraphData_.size(); ++i) {
00084     delete remoteGraphData_[i];
00085   }
00086 }
00087 
00088 //----------------------------------------------------------------------------
00089 int fei::Graph_Impl::addIndices(int row, int len, const int* indices)
00090 {
00091   if (row < 0) {
00092     return(-1);
00093   }
00094 
00095   if (row < firstLocalRow_ || row > lastLocalRow_) {
00096     int p = eqnComm_->getOwnerProc(row);
00097     remoteGraphData_[p]->addIndices(row, len, indices);
00098   }
00099   else localGraphData_->addIndices(row, len, indices);
00100 
00101   return(0);
00102 }
00103 
00104 //----------------------------------------------------------------------------
00105 int fei::Graph_Impl::addSymmetricIndices(int numIndices, int* indices,
00106           bool diagonal)
00107 {
00108   if (diagonal) {
00109     addDiagonals(numIndices, indices);
00110     return(0);
00111   }
00112 
00113   bool all_local = true;
00114   if (numProcs_ > 1) {
00115     for(int i=0; i<numIndices; ++i) {
00116       if (indices[i] < 0) {
00117   return(-1);
00118       }
00119 
00120       bool local = true;
00121       if (indices[i] < firstLocalRow_ || indices[i] > lastLocalRow_) {
00122     all_local = false;
00123     local = false;
00124       }
00125 
00126       if (!local) {
00127         int p = eqnComm_->getOwnerProc(indices[i]);
00128   remoteGraphData_[p]->addIndices(indices[i], numIndices, indices);
00129       }
00130     }
00131   }
00132 
00133   if (all_local) {
00134     localGraphData_->addIndices(numIndices, indices,
00135           numIndices, indices);
00136   }
00137   else {
00138     for(int i=0; i<numIndices; ++i) {
00139       if (indices[i] >= firstLocalRow_ && indices[i] <= lastLocalRow_) {
00140     localGraphData_->addIndices(indices[i], numIndices, indices);
00141       }
00142     }
00143   }
00144 
00145   return(0);
00146 }
00147 
00148 //----------------------------------------------------------------------------
00149 void fei::Graph_Impl::addDiagonals(int numIndices, int* indices)
00150 {
00151   bool all_local = true;
00152   int i;
00153   if (numProcs_ > 1) {
00154     for(i=0; i<numIndices; ++i) {
00155       int ind = indices[i];
00156       if (ind < 0) {
00157         throw std::runtime_error("fei::Graph_Impl::addDiagonals given negative index");
00158       }
00159 
00160       bool local = true;
00161       if (ind < firstLocalRow_ || ind > lastLocalRow_) {
00162     all_local = false;
00163     local = false;
00164       }
00165 
00166       if (!local) {
00167         int p=eqnComm_->getOwnerProc(ind);
00168   remoteGraphData_[p]->addIndices(ind, 1, &ind);
00169       }
00170     }
00171   }
00172 
00173   if (all_local) {
00174     localGraphData_->addDiagonals(numIndices, indices);
00175   }
00176   else {
00177     for(i=0; i<numIndices; ++i) {
00178       int ind = indices[i];
00179       if (ind >= firstLocalRow_ && ind <= lastLocalRow_) {
00180     localGraphData_->addIndices(ind, 1, &ind);
00181       }
00182     }
00183   }
00184 }
00185 
00186 //----------------------------------------------------------------------------
00187 int fei::Graph_Impl::writeLocalGraph(FEI_OSTREAM& os, bool debug,
00188             bool prefixLinesWithPoundSign)
00189 {
00190   if (localGraphData_ == NULL) {
00191     if (prefixLinesWithPoundSign) {
00192       os << "# fei::Graph_Impl::writeLocalGraph numRows: 0"<<FEI_ENDL;
00193     }
00194     else {
00195       os << "0" << FEI_ENDL;
00196     }
00197     return(0);
00198   }
00199 
00200   if (prefixLinesWithPoundSign) {
00201     os << "# fei::Graph_Impl::writeLocalGraph numRows: ";
00202   }
00203   os << localGraphData_->getMap().size() <<FEI_ENDL;
00204 
00205   if (prefixLinesWithPoundSign) {
00206     writeToStream(*localGraphData_, os, "# ");
00207   }
00208   else {
00209     writeToStream(*localGraphData_, os, NULL);
00210   }
00211 
00212   return( 0 );
00213 }
00214 
00215 //----------------------------------------------------------------------------
00216 int fei::Graph_Impl::writeRemoteGraph(FEI_OSTREAM& os)
00217 {
00218   for(unsigned p=0; p<remoteGraphData_.size(); ++p) {
00219     writeToStream(*remoteGraphData_[p], os, "# ");
00220   }
00221 
00222   return( 0 );
00223 }
00224 
00225 //----------------------------------------------------------------------------
00226 int fei::Graph_Impl::gatherFromOverlap()
00227 {
00228   if (numProcs_ == 1) return(0);
00229 
00230 #ifndef FEI_SER
00231   //this function gathers shared-but-not-owned data onto the
00232   //owning processors.
00233 
00234   //iterate the remoteGraphData_ array and create a list of processors
00235   //that we will be sending data to. (processors which own graph rows
00236   //that we share.)
00237   std::vector<int> sendProcs;
00238   for(unsigned i=0; i<remoteGraphData_.size(); ++i) {
00239     if ((int)i == localProc_) continue;
00240     if (remoteGraphData_[i] != NULL) {
00241       if (remoteGraphData_[i]->getMap().size() == 0) {
00242         continue;
00243       }
00244 
00245       sendProcs.push_back((int)i);
00246     }
00247   }
00248 
00249   //now we can find out which procs we'll be receiving from.
00250   std::vector<int> recvProcs;
00251   fei::mirrorProcs(comm_, sendProcs, recvProcs);
00252 
00253   //next we'll declare arrays to receive into.
00254   std::vector<std::vector<int> > recv_ints(recvProcs.size());
00255 
00256   //...and an array for the sizes of the recv buffers:
00257   std::vector<int> recv_sizes(recvProcs.size());
00258   std::vector<MPI_Request> mpiReqs(recvProcs.size());
00259   std::vector<MPI_Status> mpiStatuses(recvProcs.size());
00260 
00261   int tag1 = 11113;
00262 
00263   unsigned offset = 0;
00264   for(unsigned i=0; i<recvProcs.size(); ++i) {
00265     MPI_Irecv(&recv_sizes[i], 1, MPI_INT, recvProcs[i],
00266               tag1, comm_, &mpiReqs[i]);
00267   }
00268 
00269   //now we'll pack our to-be-sent data into buffers, and send the
00270   //sizes to the receiving procs:
00271   std::vector<std::vector<int> > send_ints(sendProcs.size());
00272 
00273   for(unsigned i=0; i<sendProcs.size(); ++i) {
00274     int proc = sendProcs[i];
00275 
00276     fei::packRaggedTable(*(remoteGraphData_[proc]), send_ints[i]);
00277 
00278     int isize = send_ints[i].size();
00279 
00280     MPI_Send(&isize, 1, MPI_INT, proc, tag1, comm_);
00281   }
00282 
00283   if (mpiReqs.size() > 0) {
00284     MPI_Waitall(mpiReqs.size(), &mpiReqs[0], &mpiStatuses[0]);
00285   }
00286 
00287   //now resize our recv buffers, and post the recvs.
00288   for(size_t i=0; i<recvProcs.size(); ++i) {
00289     int intsize = recv_sizes[i];
00290 
00291     recv_ints[i].resize(intsize);
00292 
00293     MPI_Irecv(&(recv_ints[i][0]), intsize, MPI_INT, recvProcs[i],
00294               tag1, comm_, &mpiReqs[i]);
00295   }
00296 
00297   //now send our packed buffers.
00298   for(size_t i=0; i<sendProcs.size(); ++i) {
00299     int proc = sendProcs[i];
00300 
00301     MPI_Send(&(send_ints[i][0]), send_ints[i].size(), MPI_INT,
00302              proc, tag1, comm_);
00303   }
00304 
00305   if (mpiReqs.size() > 0) {
00306     MPI_Waitall(mpiReqs.size(), &mpiReqs[0], &mpiStatuses[0]);
00307   }
00308 
00309   for(unsigned i=0; i<recvProcs.size(); ++i) {
00310     std::vector<int> recvdata = recv_ints[i];
00311     int numRows = recvdata[0];
00312     int* rowNumbers = &recvdata[1];
00313     int* rowLengths = rowNumbers+numRows;
00314     int* packedCols = rowLengths+numRows;
00315     offset = 0;
00316     for(int r=0; r<numRows; ++r) {
00317       addIndices(rowNumbers[r], rowLengths[r], &packedCols[offset]);
00318       offset += rowLengths[r];
00319     }
00320   }
00321 
00322 #endif
00323   return(0);
00324 }
00325 
00326 //----------------------------------------------------------------------------
00327 int fei::Graph_Impl::getLocalRowLength(int row)
00328 {
00329   table_row_type* colIndices = localGraphData_->getRow(row);
00330   if (colIndices == NULL) {
00331     return(-1);
00332   }
00333 
00334   return(colIndices->size());
00335 }
00336 
00337 //----------------------------------------------------------------------------
00338 int fei::Graph_Impl::getNumLocalRows()
00339 {
00340   int numLocalRows = localGraphData_->getMap().size();
00341 
00342   return(numLocalRows);
00343 }
00344 
00345 //----------------------------------------------------------------------------
00346 int fei::Graph_Impl::getNumLocalNonzeros() const
00347 {
00348   int numNonzeros = fei::countNonzeros(*localGraphData_);
00349 
00350   return(numNonzeros);
00351 }
00352 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends