Epetra Package Browser (Single Doxygen Collection) Development
Epetra_MpiDistributor.cpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 // 
00004 //               Epetra: Linear Algebra Services Package 
00005 //                 Copyright 2011 Sandia Corporation
00006 // 
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the 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 Michael A. Heroux (maherou@sandia.gov) 
00038 // 
00039 // ************************************************************************
00040 //@HEADER
00041 
00042 #include "Epetra_MpiDistributor.h"
00043 #include "Epetra_MpiComm.h"
00044 
00045 
00046 //==============================================================================
00047 Epetra_MpiDistributor::Epetra_MpiDistributor(const Epetra_MpiComm & Comm): 
00048   Epetra_Object("Epetra::MpiDistributor"),
00049   lengths_to_(0),
00050   procs_to_(0),
00051   indices_to_(0),
00052   size_indices_to_(0),
00053   lengths_from_(0),
00054   procs_from_(0),
00055   indices_from_(0),
00056   size_indices_from_(0),
00057   resized_(false),
00058   sizes_(0),
00059   sizes_to_(0),
00060   starts_to_(0),
00061   starts_to_ptr_(0),
00062   indices_to_ptr_(0),
00063   sizes_from_(0),
00064   starts_from_(0),
00065   starts_from_ptr_(0),
00066   indices_from_ptr_(0),
00067   nrecvs_(0),
00068   nsends_(0),
00069   nexports_(0),
00070   self_msg_(0),
00071   max_send_length_(0),
00072   total_recv_length_(0),
00073   tag_(Comm.GetMpiTag()),
00074   epComm_(&Comm),
00075   comm_(Comm.GetMpiComm()),
00076   request_(0),
00077   status_(0),
00078   no_delete_(false),
00079   send_array_(0),
00080   send_array_size_(0),
00081   comm_plan_reverse_(0)
00082 {
00083 }
00084 
00085 //==============================================================================
00086 Epetra_MpiDistributor::Epetra_MpiDistributor(const Epetra_MpiDistributor & Distributor):
00087   Epetra_Object("Epetra::MpiDistributor"),
00088   lengths_to_(0),
00089   procs_to_(0),
00090   indices_to_(0),
00091   size_indices_to_(Distributor.size_indices_to_),
00092   lengths_from_(0),
00093   procs_from_(0),
00094   indices_from_(0),
00095   size_indices_from_(Distributor.size_indices_from_),
00096   resized_(false),
00097   sizes_(0),
00098   sizes_to_(0),
00099   starts_to_(0),
00100   starts_to_ptr_(0),
00101   indices_to_ptr_(0),
00102   sizes_from_(0),
00103   starts_from_(0),
00104   starts_from_ptr_(0),
00105   indices_from_ptr_(0),
00106   nrecvs_(Distributor.nrecvs_),
00107   nsends_(Distributor.nsends_),
00108   nexports_(Distributor.nexports_),
00109   self_msg_(Distributor.self_msg_),
00110   max_send_length_(Distributor.max_send_length_),
00111   total_recv_length_(Distributor.total_recv_length_),
00112   tag_(Distributor.tag_),
00113   epComm_(Distributor.epComm_),
00114   comm_(Distributor.comm_),
00115   request_(0),
00116   status_(0),
00117   no_delete_(Distributor.no_delete_),
00118   send_array_(0),
00119   send_array_size_(0),
00120   comm_plan_reverse_(0)
00121 {
00122   int i;
00123   if (nsends_>0) {
00124     lengths_to_ = new int[nsends_];
00125     procs_to_ = new int[nsends_];
00126     for (i=0; i<nsends_; i++) {
00127       lengths_to_[i] = Distributor.lengths_to_[i];
00128       procs_to_[i] = Distributor.procs_to_[i];
00129     }
00130   }
00131   if (size_indices_to_>0) {
00132     indices_to_ = new int[size_indices_to_];
00133     for (i=0; i<size_indices_to_; i++) {
00134       indices_to_[i] = Distributor.indices_to_[i];
00135     }
00136   }
00137 
00138   if (nrecvs_>0) {
00139     lengths_from_ = new int[nrecvs_];
00140     procs_from_ = new int[nrecvs_];
00141     request_ = new MPI_Request[ nrecvs_ ];
00142     status_ = new MPI_Status[ nrecvs_ ];
00143     for (i=0; i<nrecvs_; i++) {
00144       lengths_from_[i] = Distributor.lengths_from_[i];
00145       procs_from_[i] = Distributor.procs_from_[i];
00146     }
00147   }
00148   if (size_indices_from_>0) {
00149     indices_from_ = new int[size_indices_from_];
00150     for (i=0; i<size_indices_from_; i++) {
00151       indices_from_[i] = Distributor.indices_from_[i];
00152     }
00153   }
00154 }
00155 
00156 //==============================================================================
00157 Epetra_MpiDistributor::~Epetra_MpiDistributor()
00158 {
00159   if( !no_delete_ )
00160   {
00161     if( lengths_to_ != 0 ) delete [] lengths_to_;
00162     if( procs_to_ != 0 ) delete [] procs_to_;
00163     if( indices_to_ != 0 ) delete [] indices_to_;
00164     if( lengths_from_ != 0 ) delete [] lengths_from_;
00165     if( procs_from_ != 0 ) delete [] procs_from_;
00166     if( indices_from_ != 0 ) delete [] indices_from_;
00167     if( starts_to_ != 0 ) delete [] starts_to_;
00168     if( starts_from_ != 0 ) delete [] starts_from_;
00169   }
00170 
00171   if( sizes_ != 0 ) delete [] sizes_;
00172   if( sizes_to_ != 0 ) delete [] sizes_to_;
00173   if( sizes_from_ != 0 ) delete [] sizes_from_;
00174   if( starts_to_ptr_ != 0 ) delete [] starts_to_ptr_;
00175   if( starts_from_ptr_ != 0 ) delete [] starts_from_ptr_;
00176   if( indices_to_ptr_ != 0 ) delete [] indices_to_ptr_;
00177   if( indices_from_ptr_ != 0 ) delete [] indices_from_ptr_;
00178 
00179   if( request_ != 0 ) delete [] request_;
00180   if( status_ != 0 ) delete [] status_;
00181 
00182   if( send_array_size_ != 0 ) { delete [] send_array_; send_array_size_ = 0; }
00183 
00184   if( comm_plan_reverse_ != 0 ) delete comm_plan_reverse_;
00185 }
00186 
00187 
00188 //==============================================================================
00189 int Epetra_MpiDistributor::CreateFromSends( const int & NumExportIDs,
00190                                             const int * ExportPIDs,
00191                                             bool Deterministic,
00192                                             int & NumRemoteIDs )
00193 {
00194   (void)Deterministic; // Prevent compiler warnings for unused argument.
00195   nexports_ = NumExportIDs;
00196 
00197   int i;
00198 
00199   int my_proc;
00200   MPI_Comm_rank( comm_, &my_proc );
00201 
00202   int nprocs;
00203   MPI_Comm_size( comm_, &nprocs );
00204 
00205   // Check to see if items are grouped by processor w/o gaps
00206   // If so, indices_to -> 0
00207 
00208   // Setup data structures for quick traversal of arrays
00209   int * starts = new int[ nprocs + 1 ];
00210   for( i = 0; i < nprocs; i++ )
00211     starts[i] = 0;
00212 
00213   int nactive = 0;
00214   bool no_send_buff = true;
00215   int numDeadIndices = 0; // In some cases the GIDs will not be owned by any processors and the PID will be -1
00216 
00217   for( i = 0; i < NumExportIDs; i++ )
00218   {
00219     if( no_send_buff && i && (ExportPIDs[i] < ExportPIDs[i-1]) )
00220       no_send_buff = false;
00221     if( ExportPIDs[i] >= 0 )
00222     {
00223       ++starts[ ExportPIDs[i] ];
00224       ++nactive;
00225     }
00226     else numDeadIndices++; // Increase the number of dead indices.  Used below to leave these out of the analysis
00227   }
00228 
00229   self_msg_ = ( starts[my_proc] != 0 ) ? 1 : 0;
00230 
00231   nsends_ = 0;
00232 
00233   if( no_send_buff ) //grouped by processor, no send buffer or indices_to_ needed
00234   {
00235     for( i = 0; i < nprocs; ++i )
00236       if( starts[i] ) ++nsends_;
00237 
00238     if( nsends_ )
00239     {
00240       procs_to_ = new int[nsends_];
00241       starts_to_ = new int[nsends_];
00242       lengths_to_ = new int[nsends_];
00243     }
00244 
00245     int index = numDeadIndices;  // Leave off the dead indices (PID = -1)
00246     int proc;
00247     for( i = 0; i < nsends_; ++i )
00248     {
00249       starts_to_[i] = index;
00250       proc = ExportPIDs[index];
00251       procs_to_[i] = proc;
00252       index += starts[proc];
00253     }
00254 
00255     if( nsends_ )
00256       Sort_ints_( procs_to_, starts_to_, nsends_ );
00257 
00258     max_send_length_ = 0;
00259 
00260     for( i = 0; i < nsends_; ++i )
00261     {
00262       proc = procs_to_[i];
00263       lengths_to_[i] = starts[proc];
00264       if( (proc != my_proc) && (lengths_to_[i] > max_send_length_) )
00265         max_send_length_ = lengths_to_[i];
00266     }
00267   }
00268   else //not grouped by processor, need send buffer and indices_to_
00269   {
00270     if( starts[0] != 0 ) nsends_ = 1;
00271 
00272     for( i = 1; i < nprocs; i++ )
00273     {
00274       if( starts[i] != 0 ) ++nsends_;
00275       starts[i] += starts[i-1];
00276     }
00277 
00278     for( i = nprocs-1; i != 0; i-- )
00279       starts[i] = starts[i-1];
00280 
00281     starts[0] = 0;
00282 
00283     if (nactive>0) {
00284       indices_to_ = new int[ nactive ];
00285       size_indices_to_ = nactive;
00286     }
00287 
00288     for( i = 0; i < NumExportIDs; i++ )
00289     if( ExportPIDs[i] >= 0 )
00290     {
00291       indices_to_[ starts[ ExportPIDs[i] ] ] = i;
00292       ++starts[ ExportPIDs[i] ];
00293     }
00294 
00295     //Reconstuct starts array to index into indices_to.
00296 
00297     for( i = nprocs-1; i != 0; i-- )
00298       starts[i] = starts[i-1];
00299     starts[0] = 0;
00300     starts[nprocs] = nactive;
00301 
00302     if (nsends_>0) {
00303       lengths_to_ = new int[ nsends_ ];
00304       procs_to_ = new int[ nsends_ ];
00305       starts_to_ = new int[ nsends_ ];
00306     }
00307 
00308     int j = 0;
00309     max_send_length_ = 0;
00310 
00311     for( i = 0; i < nprocs; i++ )
00312       if( starts[i+1] != starts[i] )
00313       {
00314         lengths_to_[j] = starts[i+1] - starts[i];
00315         starts_to_[j] = starts[i];
00316         if( ( i != my_proc ) && ( lengths_to_[j] > max_send_length_ ) )
00317           max_send_length_ = lengths_to_[j];
00318         procs_to_[j] = i;
00319         j++;
00320       }
00321   }
00322     
00323   delete [] starts;
00324 
00325   nsends_ -= self_msg_;
00326 
00327   //Invert map to see what msgs are received and what length
00328   EPETRA_CHK_ERR( ComputeRecvs_( my_proc, nprocs ) );
00329 
00330   if (nrecvs_>0) {
00331     if( !request_ ) {
00332       request_ = new MPI_Request[ nrecvs_ ];
00333       status_ = new MPI_Status[ nrecvs_ ];
00334     }
00335   }
00336 
00337   NumRemoteIDs = total_recv_length_;
00338 
00339   return 0;
00340 }
00341 
00342 //==============================================================================
00343 int Epetra_MpiDistributor::CreateFromRecvs( const int & NumRemoteIDs,
00344            const int * RemoteGIDs,
00345                  const int * RemotePIDs,
00346            bool Deterministic,
00347                  int & NumExportIDs,
00348            int *& ExportGIDs,
00349            int *& ExportPIDs )
00350 {
00351   int my_proc;
00352   MPI_Comm_rank( comm_, &my_proc );
00353 
00354   int nprocs;
00355   MPI_Comm_size( comm_, &nprocs );
00356 
00357   EPETRA_CHK_ERR( ComputeSends_( NumRemoteIDs, RemoteGIDs, RemotePIDs, NumExportIDs,
00358          ExportGIDs, ExportPIDs, my_proc) );
00359 
00360   int testNumRemoteIDs;
00361   EPETRA_CHK_ERR( CreateFromSends( NumExportIDs, ExportPIDs,
00362            Deterministic, testNumRemoteIDs ) );
00363 
00364   return(0);
00365 }
00366 
00367 //==============================================================================
00368 //---------------------------------------------------------------------------
00369 //ComputeRecvs Method
00370 //---------------------------------------------------------------------------
00371 int Epetra_MpiDistributor::ComputeRecvs_( int my_proc, 
00372               int nprocs )
00373 {
00374   int * msg_count = new int[ nprocs ];
00375   int * counts = new int[ nprocs ];
00376 
00377   int i;
00378   MPI_Status status;
00379 
00380   for( i = 0; i < nprocs; i++ )
00381   {
00382     msg_count[i] = 0;
00383     counts[i] = 1;
00384   }
00385 
00386   for( i = 0; i < nsends_+self_msg_; i++ )
00387     msg_count[ procs_to_[i] ] = 1;
00388 
00389 #if defined(REDUCE_SCATTER_BUG)
00390 // the bug is found in mpich on linux platforms
00391   MPI_Reduce(msg_count, counts, nprocs, MPI_INT, MPI_SUM, 0, comm_);
00392   MPI_Scatter(counts, 1, MPI_INT, &nrecvs_, 1, MPI_INT, 0, comm_);
00393 #else
00394   MPI_Reduce_scatter( msg_count, &nrecvs_, counts, MPI_INT, MPI_SUM, comm_ );
00395 #endif
00396 
00397   delete [] msg_count;
00398   delete [] counts;
00399 
00400   if (nrecvs_>0) {
00401     lengths_from_ = new int[nrecvs_];
00402     procs_from_ = new int[nrecvs_];
00403     for(i=0; i<nrecvs_; ++i) {
00404       lengths_from_[i] = 0;
00405       procs_from_[i] = 0;
00406     }
00407   }
00408 
00409 #ifndef NEW_COMM_PATTERN
00410   for( i = 0; i < (nsends_+self_msg_); i++ )
00411     if( procs_to_[i] != my_proc ) {
00412       MPI_Send( &(lengths_to_[i]), 1, MPI_INT, procs_to_[i], tag_, comm_ );
00413     }
00414     else
00415     {
00416       //set self_msg_ to end block of recv arrays
00417       lengths_from_[nrecvs_-1] = lengths_to_[i];
00418       procs_from_[nrecvs_-1] = my_proc;
00419     }
00420 
00421   for( i = 0; i < (nrecvs_-self_msg_); i++ )
00422   {
00423     MPI_Recv( &(lengths_from_[i]), 1, MPI_INT, MPI_ANY_SOURCE, tag_, comm_, &status );
00424     procs_from_[i] = status.MPI_SOURCE;
00425   }
00426 
00427   MPI_Barrier( comm_ );
00428 #else
00429   if (nrecvs_>0) {
00430     if( !request_ ) {
00431       request_ = new MPI_Request[nrecvs_-self_msg_];
00432       status_ = new MPI_Status[nrecvs_-self_msg_];
00433     }
00434   }
00435 
00436   for( i = 0; i < (nrecvs_-self_msg_); i++ )
00437     MPI_Irecv( &(lengths_from_[i]), 1, MPI_INT, MPI_ANY_SOURCE, tag_, comm_, &(request_[i]) );
00438 
00439   MPI_Barrier( comm_ );
00440 
00441   for( i = 0; i < (nsends_+self_msg_); i++ )
00442     if( procs_to_[i] != my_proc ) {
00443       MPI_Rsend( &(lengths_to_[i]), 1, MPI_INT, procs_to_[i], tag_, comm_ );
00444     }
00445     else
00446     {
00447       //set self_msg_ to end block of recv arrays
00448       lengths_from_[nrecvs_-1] = lengths_to_[i];
00449       procs_from_[nrecvs_-1] = my_proc;
00450     }
00451 
00452   if( (nrecvs_-self_msg_) > 0 ) MPI_Waitall( (nrecvs_-self_msg_), request_, status_ );
00453 
00454   for( i = 0; i < (nrecvs_-self_msg_); i++ )
00455     procs_from_[i] = status_[i].MPI_SOURCE;
00456 #endif
00457 
00458   Sort_ints_( procs_from_, lengths_from_, nrecvs_ );
00459 
00460   // Compute indices_from_
00461   // Seems to break some rvs communication
00462 /* Not necessary since rvs communication is always blocked
00463   size_indices_from_ = 0;
00464   if( nrecvs_ > 0 )
00465   {
00466     for( i = 0; i < nrecvs_; i++ )  size_indices_from_ += lengths_from_[i];
00467     indices_from_ = new int[ size_indices_from_ ];
00468 
00469     for (i=0; i<size_indices_from_; i++) indices_from_[i] = i;
00470   }
00471 */
00472 
00473   if (nrecvs_>0) starts_from_ = new int[nrecvs_];
00474   int j = 0;
00475   for( i=0; i<nrecvs_; ++i )
00476   {
00477     starts_from_[i] = j;
00478     j += lengths_from_[i];
00479   }
00480 
00481   total_recv_length_ = 0;
00482   for( i = 0; i < nrecvs_; i++ )
00483     total_recv_length_ += lengths_from_[i];
00484 
00485   nrecvs_ -= self_msg_;
00486 
00487   MPI_Barrier( comm_ );
00488   
00489   return false;
00490 }
00491 
00492 //==============================================================================
00493 //---------------------------------------------------------------------------
00494 //ComputeSends Method
00495 //---------------------------------------------------------------------------
00496 int Epetra_MpiDistributor::ComputeSends_( int num_imports,
00497         const int *& import_ids,
00498         const int *& import_procs,
00499         int & num_exports,
00500         int *& export_ids,
00501         int *& export_procs,
00502         int my_proc ) {
00503  
00504   Epetra_MpiDistributor tmp_plan(*epComm_);
00505   int i;
00506 
00507   int * proc_list = 0;
00508   int * import_objs = 0;
00509   char * c_export_objs = 0;
00510 
00511   if( num_imports > 0 )
00512   {
00513     proc_list = new int[ num_imports ];
00514     import_objs = new int[ 2 * num_imports ];
00515 
00516     for( i = 0; i < num_imports; i++ )
00517     {
00518       proc_list[i] = import_procs[i];
00519 
00520       import_objs[2*i] = import_ids[i];
00521       import_objs[2*i+1] = my_proc;
00522     }
00523   }
00524 
00525   EPETRA_CHK_ERR(tmp_plan.CreateFromSends( num_imports, proc_list,
00526              true, num_exports) );
00527   if( num_exports > 0 )
00528   {
00529     //export_objs = new int[ 2 * num_exports ];
00530     export_ids = new int[ num_exports ];
00531     export_procs = new int[ num_exports ];
00532   }
00533   else
00534   {
00535     export_ids = 0;
00536     export_procs = 0;
00537   }
00538 
00539   int len_c_export_objs = 0;
00540   EPETRA_CHK_ERR( tmp_plan.Do(reinterpret_cast<char *> (import_objs),
00541             2 * (int)sizeof( int ), 
00542             len_c_export_objs,
00543             c_export_objs) );
00544   int * export_objs = reinterpret_cast<int *>(c_export_objs);
00545 
00546   for( i = 0; i < num_exports; i++ ) {
00547     export_ids[i] = export_objs[2*i];
00548     export_procs[i] = export_objs[2*i+1];
00549   }
00550 
00551   if( proc_list != 0 ) delete [] proc_list;
00552   if( import_objs != 0 ) delete [] import_objs;
00553   if( len_c_export_objs != 0 ) delete [] c_export_objs;
00554 
00555   return(0);
00556 
00557 }
00558 
00559 //==============================================================================
00560 int Epetra_MpiDistributor::Do( char * export_objs,
00561                                int obj_size,
00562                                int & len_import_objs,
00563                                char *& import_objs )
00564 {
00565   EPETRA_CHK_ERR( DoPosts(export_objs, obj_size, len_import_objs, import_objs) );
00566   EPETRA_CHK_ERR( DoWaits() );
00567   return(0);
00568 }
00569 
00570 //==============================================================================
00571 int Epetra_MpiDistributor::DoReverse( char * export_objs,
00572                                       int obj_size,
00573                                       int & len_import_objs,
00574                                       char *& import_objs )
00575 {
00576   EPETRA_CHK_ERR( DoReversePosts(export_objs, obj_size,
00577          len_import_objs, import_objs) );
00578   EPETRA_CHK_ERR( DoReverseWaits() );
00579   return(0);
00580 }
00581 
00582 //==============================================================================
00583 int Epetra_MpiDistributor::DoPosts( char * export_objs,
00584             int obj_size,
00585                                     int & len_import_objs,
00586             char *& import_objs )
00587 {
00588   int i, j, k;
00589 
00590   int my_proc = 0;
00591   int self_recv_address = 0;
00592 
00593   MPI_Comm_rank( comm_, &my_proc );
00594 
00595   if( len_import_objs < (total_recv_length_*obj_size) )
00596   {
00597     if( import_objs!=0 ) {delete [] import_objs; import_objs = 0;}
00598     len_import_objs = total_recv_length_*obj_size;
00599     if (len_import_objs>0) import_objs = new char[len_import_objs];
00600     for( i=0; i<len_import_objs; ++i ) import_objs[i]=0;
00601   }
00602 
00603   k = 0;
00604 
00605   j = 0;
00606   for( i = 0; i < (nrecvs_+self_msg_); i++ )
00607   {
00608     if( procs_from_[i] != my_proc )
00609     {
00610       MPI_Irecv( &(import_objs[j]),
00611                  lengths_from_[i] * obj_size,
00612                  MPI_CHAR, procs_from_[i],
00613                  tag_, comm_,
00614                  &(request_[k]) );
00615       k++;
00616     }
00617     else
00618       self_recv_address = j;
00619 
00620     j += lengths_from_[i] * obj_size;
00621   }
00622 
00623 #ifndef EPETRA_NO_READY_SEND_IN_DO_POSTS
00624   // NOTE (mfh 19 Mar 2012):
00625   //
00626   // The ready-sends below require that each ready-send's matching
00627   // receive (see above) has already been posted.  We ensure this with
00628   // a barrier.  (Otherwise, some process that doesn't need to post
00629   // receives might post its ready-send before the receiving process
00630   // gets to post its receive.)  If you want to remove the barrier,
00631   // you'll have to replace the ready-sends below with standard sends
00632   // or Isends.
00633   MPI_Barrier( comm_ );
00634 #endif // EPETRA_NO_READY_SEND_IN_DO_POSTS
00635 
00636   //setup scan through procs_to list starting w/ higher numbered procs 
00637   //Should help balance msg traffic
00638   int nblocks = nsends_ + self_msg_; 
00639   int proc_index = 0; 
00640   while( proc_index < nblocks && procs_to_[proc_index] < my_proc )
00641     ++proc_index;                    
00642   if( proc_index == nblocks ) proc_index = 0;
00643    
00644   int self_num = 0, self_index = 0;
00645   int p;
00646 
00647   if( !indices_to_ ) //data already blocked by processor
00648   {
00649     for( i = 0; i < nblocks; ++i )
00650     {
00651       p = i + proc_index;
00652       if( p > (nblocks-1) ) p -= nblocks;
00653 
00654       if( procs_to_[p] != my_proc ) {
00655 #ifndef EPETRA_NO_READY_SEND_IN_DO_POSTS
00656         MPI_Rsend( &export_objs[starts_to_[p]*obj_size],
00657                    lengths_to_[p]*obj_size,
00658                    MPI_CHAR,
00659                    procs_to_[p],
00660                    tag_,
00661                    comm_ );
00662 #else
00663         MPI_Send( &export_objs[starts_to_[p]*obj_size],
00664                    lengths_to_[p]*obj_size,
00665                    MPI_CHAR,
00666                    procs_to_[p],
00667                    tag_,
00668                    comm_ );
00669 #endif // EPETRA_NO_READY_SEND_IN_DO_POSTS
00670       }
00671       else {
00672        self_num = p;
00673       }
00674     }
00675 
00676     if( self_msg_ )
00677       memcpy( &import_objs[self_recv_address],
00678               &export_objs[starts_to_[self_num]*obj_size],
00679               lengths_to_[self_num]*obj_size );
00680   }
00681   else //data not blocked by proc, use send buffer
00682   {
00683     if( send_array_size_ < (max_send_length_*obj_size) )
00684     {
00685       if( send_array_!=0 ) {delete [] send_array_; send_array_ = 0;}
00686       send_array_size_ = max_send_length_*obj_size;
00687       if (send_array_size_>0) send_array_ = new char[send_array_size_];
00688     }
00689 
00690     j = 0;
00691     for( i = 0; i < nblocks; i++ )
00692     {
00693       p = i + proc_index;
00694       if( p > (nblocks-1) ) p -= nblocks;
00695       if( procs_to_[p] != my_proc )
00696       {
00697         int offset = 0;
00698         j = starts_to_[p];
00699         for( k = 0; k < lengths_to_[p]; k++ )
00700         {
00701           memcpy( &(send_array_[offset]), 
00702                   &(export_objs[indices_to_[j]*obj_size]),
00703                   obj_size );
00704           ++j;
00705           offset += obj_size;
00706         }
00707 #ifndef EPETRA_NO_READY_SEND_IN_DO_POSTS
00708         MPI_Rsend( send_array_,
00709                    lengths_to_[p] * obj_size,
00710                    MPI_CHAR,
00711                    procs_to_[p],
00712                    tag_, comm_ );
00713 #else
00714         MPI_Send( send_array_,
00715       lengths_to_[p] * obj_size,
00716       MPI_CHAR,
00717       procs_to_[p],
00718       tag_, comm_ );
00719 #endif // EPETRA_NO_READY_SEND_IN_DO_POSTS
00720       }
00721       else
00722       {
00723         self_num = p;
00724         self_index = starts_to_[p];
00725       }
00726     }
00727 
00728     if( self_msg_ )
00729       for( k = 0; k < lengths_to_[self_num]; k++ )
00730       {
00731         memcpy( &(import_objs[self_recv_address]),
00732                 &(export_objs[indices_to_[self_index]*obj_size]),
00733                 obj_size );
00734         self_index++;
00735         self_recv_address += obj_size;
00736       }
00737   }
00738 
00739   return(0);
00740 }
00741 
00742 //==============================================================================
00743 int Epetra_MpiDistributor::DoWaits()
00744 {
00745   if( nrecvs_ > 0 ) MPI_Waitall( nrecvs_, request_, status_ );
00746 
00747   return(0);
00748 }
00749 
00750 //==============================================================================
00751 int Epetra_MpiDistributor::DoReversePosts( char * export_objs,
00752                                            int obj_size,
00753                                            int & len_import_objs,
00754                                            char *& import_objs )
00755 {
00756   assert(indices_to_==0); //Can only do reverse comm when original data
00757                           // is blocked by processor
00758 
00759   int i;
00760   int my_proc = 0;
00761 
00762   MPI_Comm_rank( comm_, &my_proc );
00763 
00764   if( comm_plan_reverse_ == 0 )
00765   {
00766     int total_send_length = 0;
00767     for( i = 0; i < nsends_+self_msg_; i++ )
00768       total_send_length += lengths_to_[i];
00769 
00770     int max_recv_length = 0;
00771     for( i = 0; i < nrecvs_; i++ )
00772       if( procs_from_[i] != my_proc )
00773         if( lengths_from_[i] > max_recv_length )
00774           max_recv_length = lengths_from_[i];
00775 
00776     comm_plan_reverse_ = new Epetra_MpiDistributor(*epComm_);
00777 
00778     comm_plan_reverse_->lengths_to_ = lengths_from_;
00779     comm_plan_reverse_->procs_to_ = procs_from_;
00780     comm_plan_reverse_->indices_to_ = indices_from_;
00781     comm_plan_reverse_->starts_to_ = starts_from_;
00782 
00783     comm_plan_reverse_->lengths_from_ = lengths_to_;
00784     comm_plan_reverse_->procs_from_ = procs_to_;
00785     comm_plan_reverse_->indices_from_ = indices_to_;
00786     comm_plan_reverse_->starts_from_ = starts_to_;
00787 
00788     comm_plan_reverse_->nsends_ = nrecvs_;
00789     comm_plan_reverse_->nrecvs_ = nsends_;
00790     comm_plan_reverse_->self_msg_ = self_msg_;
00791 
00792     comm_plan_reverse_->max_send_length_ = max_recv_length;
00793     comm_plan_reverse_->total_recv_length_ = total_send_length;
00794 
00795     comm_plan_reverse_->request_ = new MPI_Request[ comm_plan_reverse_->nrecvs_ ];
00796     comm_plan_reverse_->status_= new MPI_Status[ comm_plan_reverse_->nrecvs_ ];
00797 
00798     comm_plan_reverse_->no_delete_ = true;
00799   }
00800 
00801   int comm_flag = comm_plan_reverse_->DoPosts(export_objs, obj_size, len_import_objs, import_objs);
00802 
00803   return(comm_flag);
00804 }
00805 
00806 //==============================================================================
00807 int Epetra_MpiDistributor::DoReverseWaits()
00808 {
00809   if( comm_plan_reverse_ == 0 ) return (-1);
00810 
00811   int comm_flag = comm_plan_reverse_->DoWaits();
00812 
00813   return(comm_flag);
00814 }
00815 
00816 //==============================================================================
00817 //---------------------------------------------------------------------------
00818 //Resize Method                      (Heaphy) 
00819 //---------------------------------------------------------------------------
00820 int Epetra_MpiDistributor::Resize_( int * sizes )
00821 { 
00822   int  i, j, k;         // loop counters
00823   int sum; 
00824  
00825   //if (sizes == 0) return 0; 
00826      
00827   int my_proc; 
00828   MPI_Comm_rank (comm_, &my_proc);  
00829   int nprocs;
00830   MPI_Comm_size( comm_, &nprocs );
00831      
00832   if( resized_ )
00833   {  
00834     //test and see if we are already setup for these sizes
00835     bool match = true; 
00836     for( i = 0; i < nexports_; ++i )
00837       match = match && (sizes_[i]==sizes[i]);
00838     int matched = match?1:0;
00839     int match_count = 0;
00840     MPI_Allreduce( &matched, &match_count, 1, MPI_INT, MPI_SUM, comm_ );
00841     if( match_count == nprocs )
00842       return 0;
00843     else //reset existing sizing arrays
00844       max_send_length_ = 0;
00845   } 
00846  
00847   if( !sizes_ && nexports_ ) sizes_ = new int[nexports_]; 
00848   for (i = 0; i < nexports_; i++)
00849     sizes_[i] = sizes[i]; 
00850  
00851   if( !sizes_to_ && (nsends_+self_msg_) ) sizes_to_ = new int[nsends_+self_msg_];
00852   for (i = 0; i < (nsends_+self_msg_); ++i) 
00853     sizes_to_[i] = 0;                       
00854  
00855   if( !starts_to_ptr_ && (nsends_+self_msg_) ) starts_to_ptr_ = new int[nsends_+self_msg_];
00856 
00857   if( !indices_to_ ) //blocked sends
00858   {
00859     
00860     int * index = 0;
00861     int * sort_val = 0;
00862     if (nsends_+self_msg_>0) {
00863       index = new int[nsends_+self_msg_];
00864       sort_val = new int[nsends_+self_msg_];
00865     }
00866     for( i = 0; i < (nsends_+self_msg_); ++i )
00867     {
00868       j = starts_to_[i];
00869       for( k = 0; k < lengths_to_[i]; ++k )
00870         sizes_to_[i] += sizes[j++];
00871       if( (sizes_to_[i] > max_send_length_) && (procs_to_[i] != my_proc) )
00872         max_send_length_ = sizes_to_[i];
00873     }
00874 
00875     for( i = 0; i < (nsends_+self_msg_); ++i )
00876     {
00877       sort_val[i] = starts_to_[i];
00878       index[i] = i;
00879     }
00880 
00881     if( nsends_+self_msg_ )
00882       Sort_ints_( sort_val, index, (nsends_+self_msg_) );
00883 
00884     sum = 0;
00885     for( i = 0; i < (nsends_+self_msg_); ++i )
00886     {
00887       starts_to_ptr_[ index[i] ] = sum;
00888       sum += sizes_to_[ index[i] ];
00889     }
00890 
00891     if (index!=0) {delete [] index; index = 0;}
00892     if (sort_val!=0) {delete [] sort_val; sort_val = 0;}
00893   }
00894   else //Sends not blocked, so have to do more work
00895   {
00896     if( !indices_to_ptr_ && nexports_ ) indices_to_ptr_ = new int[nexports_]; 
00897     int * offset = 0;
00898     if( nexports_ ) offset = new int[nexports_];
00899    
00900     //Compute address for every item in send array
00901     sum = 0; 
00902     for( i = 0; i < nexports_; ++i )
00903     {  
00904       offset[i] = sum; 
00905       sum += sizes_[i];
00906     }
00907    
00908     sum = 0;
00909     max_send_length_ = 0;
00910     for( i = 0; i < (nsends_+self_msg_); ++i )
00911     {
00912       starts_to_ptr_[i] = sum;
00913       for( j = starts_to_[i]; j < (starts_to_[i]+lengths_to_[i]); ++j )
00914       {
00915         indices_to_ptr_[j] = offset[ indices_to_[j] ];
00916         sizes_to_[i] += sizes_[ indices_to_[j] ];
00917       }
00918       if( sizes_to_[i] > max_send_length_ && procs_to_[i] != my_proc )
00919         max_send_length_ = sizes_to_[i];
00920       sum += sizes_to_[i];
00921     }
00922  
00923     if (offset!=0) {delete [] offset; offset = 0;}
00924   }
00925  
00926   //  Exchange sizes routine inserted here:
00927   int self_index_to = -1;
00928   total_recv_length_ = 0;
00929   if( !sizes_from_ && (nrecvs_+self_msg_) ) sizes_from_ = new int [nrecvs_+self_msg_];
00930 
00931 #ifndef EPETRA_NEW_COMM_PATTERN
00932   for (i = 0; i < (nsends_+self_msg_); i++)
00933   {
00934     if(procs_to_[i] != my_proc)
00935       MPI_Send ((void *) &(sizes_to_[i]), 1, MPI_INT, procs_to_[i], tag_, comm_);
00936     else
00937       self_index_to = i;
00938   }
00939 
00940   MPI_Status status;
00941   for (i = 0; i < (nrecvs_+self_msg_); ++i)
00942   {
00943     sizes_from_[i] = 0;
00944     if (procs_from_[i] != my_proc)
00945       MPI_Recv((void *) &(sizes_from_[i]), 1, MPI_INT, procs_from_[i], tag_, comm_, &status);
00946     else
00947       sizes_from_[i] = sizes_to_[self_index_to];
00948     total_recv_length_ += sizes_from_[i];
00949   }
00950 #else
00951   if (nrecvs_>0 && !request_) {
00952     request_ = new MPI_Request[ nrecvs_-self_msg_ ];
00953     status_ = new MPI_Status[ nrecvs_-self_msg_ ];
00954   }
00955 
00956   for (i = 0; i < (nsends_+self_msg_); i++)
00957   {
00958     if(procs_to_[i] == my_proc)
00959       self_index_to = i;
00960   }
00961 
00962   for (i = 0; i < (nrecvs_+self_msg_); ++i)
00963   {
00964     sizes_from_[i] = 0;
00965     if (procs_from_[i] != my_proc)
00966       MPI_Irecv((void *) &(sizes_from_[i]), 1, MPI_INT, procs_from_[i], tag_, comm_, &(request_[i]));
00967     else
00968     {
00969       sizes_from_[i] = sizes_to_[self_index_to];
00970       total_recv_length_ += sizes_from_[i];
00971     }
00972   }
00973 
00974   MPI_Barrier( comm_ );
00975 
00976   for (i = 0; i < (nsends_+self_msg_); i++)
00977   {
00978     if(procs_to_[i] != my_proc)
00979       MPI_Rsend ((void *) &(sizes_to_[i]), 1, MPI_INT, procs_to_[i], tag_, comm_);
00980   }
00981 
00982   if( nrecvs_ > 0 ) MPI_Waitall( nrecvs_, request_, status_ );
00983 
00984   for (i = 0; i < (nrecvs_+self_msg_); ++i)
00985   {
00986     if (procs_from_[i] != my_proc)
00987       total_recv_length_ += sizes_from_[i];
00988   }
00989 #endif
00990   // end of exchanges sizes insert
00991  
00992   sum = 0;
00993   if( !starts_from_ptr_ ) starts_from_ptr_  = new int[nrecvs_+self_msg_]; 
00994   for (i = 0; i < (nrecvs_+self_msg_); ++i)
00995   {
00996      starts_from_ptr_[i] = sum;
00997      sum += sizes_from_[i];
00998   }
00999 
01000   resized_ = true;
01001 
01002   return 0;
01003 }
01004 
01005 //==============================================================================
01006 int Epetra_MpiDistributor::Do( char * export_objs,
01007                                int obj_size,
01008                                int *& sizes,
01009                                int & len_import_objs,
01010                                char *& import_objs )
01011 {
01012   EPETRA_CHK_ERR( DoPosts(export_objs, obj_size, sizes,
01013         len_import_objs, import_objs) );
01014   EPETRA_CHK_ERR( DoWaits() );
01015 
01016   return(0);
01017 }
01018 
01019 //==============================================================================
01020 int Epetra_MpiDistributor::DoReverse( char * export_objs,
01021                                       int obj_size,
01022                                       int *& sizes,
01023                                       int & len_import_objs,
01024                                       char *& import_objs )
01025 {
01026   EPETRA_CHK_ERR( DoReversePosts(export_objs, obj_size, sizes,
01027          len_import_objs, import_objs) );
01028   EPETRA_CHK_ERR( DoReverseWaits() );
01029 
01030   return(0);
01031 }
01032 
01033 //==============================================================================
01034 int Epetra_MpiDistributor::DoPosts( char * export_objs,
01035                                     int obj_size,
01036                                     int *& sizes,
01037                                     int & len_import_objs,
01038                                     char *& import_objs )
01039 {
01040   int ierr = Resize_(sizes);
01041   if (ierr != 0) {
01042     return(ierr);
01043   }
01044 
01045   MPI_Barrier( comm_ );
01046 
01047   int i, j, k;
01048 
01049   int my_proc = 0;
01050   int self_recv_address = 0;
01051 
01052   MPI_Comm_rank( comm_, &my_proc );
01053 
01054   if( len_import_objs < (total_recv_length_*obj_size) )
01055   {
01056     if( import_objs!=0 ) {delete [] import_objs; import_objs = 0;}
01057     len_import_objs = total_recv_length_*obj_size;
01058     if (len_import_objs>0) import_objs = new char[len_import_objs];
01059   }
01060 
01061   k = 0;
01062 
01063   for( i = 0; i < (nrecvs_+self_msg_); ++i )
01064   {
01065     if( procs_from_[i] != my_proc )
01066     {
01067       MPI_Irecv( &(import_objs[starts_from_ptr_[i] * obj_size]),
01068                  sizes_from_[i] * obj_size,
01069                  MPI_CHAR, procs_from_[i], tag_, comm_, &(request_[k]) );
01070       k++;
01071     }
01072     else
01073       self_recv_address = starts_from_ptr_[i] * obj_size;
01074   }
01075 
01076   MPI_Barrier( comm_ );
01077 
01078   //setup scan through procs_to list starting w/ higher numbered procs 
01079   //Should help balance msg traffic
01080   int nblocks = nsends_ + self_msg_; 
01081   int proc_index = 0; 
01082   while( proc_index < nblocks && procs_to_[proc_index] < my_proc )
01083     ++proc_index;                    
01084   if( proc_index == nblocks ) proc_index = 0;
01085    
01086   int self_num = 0;
01087   int p;
01088 
01089   if( !indices_to_ ) //data already blocked by processor
01090   {
01091     for( i = 0; i < nblocks; ++i )
01092     {
01093       p = i + proc_index;
01094       if( p > (nblocks-1) ) p -= nblocks;
01095 
01096       if( procs_to_[p] != my_proc )
01097         MPI_Rsend( &export_objs[starts_to_ptr_[p]*obj_size],
01098                    sizes_to_[p]*obj_size,
01099                    MPI_CHAR,
01100                    procs_to_[p],
01101                    tag_,
01102                    comm_ );
01103       else
01104         self_num = p;
01105     }
01106 
01107     if( self_msg_ )
01108       memcpy( &import_objs[self_recv_address],
01109               &export_objs[starts_to_ptr_[self_num]*obj_size],
01110               sizes_to_[self_num]*obj_size );
01111   }
01112   else //data not blocked by proc, need to copy to buffer
01113   {
01114     if( send_array_size_ && send_array_size_ < (max_send_length_*obj_size) )
01115     {
01116       if (send_array_!=0) {delete [] send_array_; send_array_ = 0;}
01117       send_array_ = 0;
01118       send_array_size_ = 0;
01119     }
01120     if( !send_array_size_ )
01121     {
01122       send_array_size_ = max_send_length_*obj_size;
01123       if (send_array_size_>0) send_array_ = new char[send_array_size_];
01124     }
01125 
01126     for( i=0; i<nblocks; ++i )
01127     {
01128       p = i + proc_index;
01129       if( p > (nblocks-1) ) p -= nblocks;
01130 
01131       if( procs_to_[p] != my_proc )
01132       {
01133         int offset = 0;
01134         j = starts_to_[p];
01135         for( k=0; k<lengths_to_[p]; ++k )
01136         {
01137           memcpy( &send_array_[offset],
01138                   &export_objs[indices_to_ptr_[j]*obj_size],
01139                   sizes_[indices_to_[j]]*obj_size );
01140           offset += sizes_[indices_to_[j]]*obj_size;
01141           ++j;
01142         }
01143         MPI_Rsend( send_array_, sizes_to_[p]*obj_size,
01144                    MPI_CHAR, procs_to_[p], tag_, comm_ );
01145       }
01146       else
01147         self_num = p;
01148     }
01149 
01150     if( self_msg_ )
01151     {
01152       int jj;
01153       j = starts_to_[self_num];
01154       for( k=0; k<lengths_to_[self_num]; ++k )
01155       {
01156         jj = indices_to_ptr_[j];
01157         memcpy( &import_objs[self_recv_address],
01158                 &export_objs[jj*obj_size],
01159                 sizes_[indices_to_[j]*obj_size] );
01160         self_recv_address += (obj_size*sizes_[indices_to_[j]]);
01161         ++jj;
01162       }
01163     }
01164   }
01165 
01166   return(0);
01167 }
01168 
01169 //==============================================================================
01170 int Epetra_MpiDistributor::DoReversePosts( char * export_objs,
01171                    int obj_size,
01172                                            int *& sizes,
01173                                            int & len_import_objs,
01174                    char *& import_objs )
01175 {
01176   assert(indices_to_==0); //Can only do reverse comm when original data
01177                           // is blocked by processor
01178 
01179   int i;
01180   int my_proc = 0;
01181 
01182   MPI_Comm_rank( comm_, &my_proc );
01183 
01184   if( comm_plan_reverse_ == 0 )
01185   {
01186     int total_send_length = 0;
01187     for( i = 0; i < nsends_+self_msg_; i++ )
01188       total_send_length += lengths_to_[i];
01189 
01190     int max_recv_length = 0;
01191     for( i = 0; i < nrecvs_; i++ )
01192       if( procs_from_[i] != my_proc )
01193         if( lengths_from_[i] > max_recv_length )
01194           max_recv_length = lengths_from_[i];
01195 
01196     comm_plan_reverse_ = new Epetra_MpiDistributor(*epComm_);
01197 
01198     comm_plan_reverse_->lengths_to_ = lengths_from_;
01199     comm_plan_reverse_->procs_to_ = procs_from_;
01200     comm_plan_reverse_->indices_to_ = indices_from_;
01201     comm_plan_reverse_->starts_to_ = starts_from_;
01202 
01203     comm_plan_reverse_->lengths_from_ = lengths_to_;
01204     comm_plan_reverse_->procs_from_ = procs_to_;
01205     comm_plan_reverse_->indices_from_ = indices_to_;
01206     comm_plan_reverse_->starts_from_ = starts_to_;
01207 
01208     comm_plan_reverse_->nsends_ = nrecvs_;
01209     comm_plan_reverse_->nrecvs_ = nsends_;
01210     comm_plan_reverse_->self_msg_ = self_msg_;
01211 
01212     comm_plan_reverse_->max_send_length_ = max_recv_length;
01213     comm_plan_reverse_->total_recv_length_ = total_send_length;
01214 
01215     comm_plan_reverse_->request_ = new MPI_Request[ comm_plan_reverse_->nrecvs_ ];
01216     comm_plan_reverse_->status_= new MPI_Status[ comm_plan_reverse_->nrecvs_ ];
01217 
01218     comm_plan_reverse_->no_delete_ = true;
01219   }
01220 
01221   int comm_flag = comm_plan_reverse_->DoPosts(export_objs, obj_size, sizes, len_import_objs, import_objs);
01222 
01223   return(comm_flag);
01224 }
01225 
01226 //==============================================================================
01227 void Epetra_MpiDistributor::Print( ostream & os) const
01228 {
01229   using std::endl;
01230 
01231   int myRank = 0, numProcs = 1;
01232   MPI_Comm_rank (comm_, &myRank);
01233   MPI_Comm_size (comm_, &numProcs);
01234 
01235   if (myRank == 0) {
01236     os << "Epetra_MpiDistributor (implements Epetra_Distributor)" << endl;
01237   }
01238   // Let each MPI process print its data.  We assume that all
01239   // processes can print to the given output stream, and execute
01240   // barriers to make it more likely that the output will be in the
01241   // right order.
01242   for (int p = 0; p < numProcs; ++p) {
01243     if (myRank == p) {
01244       os << "[Node " << p << " of " << numProcs << "]" << endl;
01245       os << " selfMessage: " << self_msg_ << endl;
01246       os << " numSends: " << nsends_ << endl;
01247 
01248       os << " imagesTo: [";
01249       for (int i = 0; i < nsends_; ++i) {
01250   os << procs_to_[i];
01251   if (i < nsends_ - 1) {
01252     os << " ";
01253   }
01254       }
01255       os << "]" << endl;
01256 
01257       os << " lengthsTo: [";
01258       for (int i = 0; i < nsends_; ++i) {
01259   os << lengths_to_[i];
01260   if (i < nsends_ - 1) {
01261     os << " ";
01262   }
01263       }
01264       os << "]" << endl;
01265 
01266       os << " maxSendLength: " << max_send_length_ << endl;
01267 
01268       os << " startsTo: ";
01269       if (starts_to_ == NULL) {
01270   os << "(NULL)" << endl;
01271       } else {
01272   os << "[";
01273   for (int i = 0; i < nsends_; ++i) {
01274     os << starts_to_[i];
01275     if (i < nsends_ - 1) {
01276       os << " ";
01277     }
01278   }
01279   os << "]" << endl;
01280       }
01281 
01282       os << " indicesTo: ";
01283       if (indices_to_ == NULL) {
01284   os << "(NULL)" << endl;
01285       } else {
01286   os << "[";
01287   int k = 0;
01288   for (int i = 0; i < nsends_; ++i) {
01289     for (int j = 0; j < lengths_to_[i]; ++j) {
01290       os << " " << indices_to_[j+k];
01291     }
01292     k += lengths_to_[i];
01293   }
01294   os << "]" << endl;
01295       }
01296 
01297       os << " numReceives: " << nrecvs_ << endl;
01298       os << " totalReceiveLength: " << total_recv_length_ << endl;
01299 
01300       os << " lengthsFrom: [";
01301       for (int i = 0; i < nrecvs_; ++i) {
01302   os << lengths_from_[i];
01303   if (i < nrecvs_ - 1) {
01304     os << " ";
01305   }
01306       }
01307       os << "]" << endl;
01308 
01309       os << " startsFrom: [";
01310       for (int i = 0; i < nrecvs_; ++i) {
01311   os << starts_from_[i];
01312   if (i < nrecvs_ - 1) {
01313     os << " ";
01314   }
01315       }
01316       os << "]" << endl;
01317 
01318       os << " imagesFrom: [";
01319       for (int i = 0; i < nrecvs_; ++i) {
01320   os << procs_from_[i];
01321   if (i < nrecvs_ - 1) {
01322     os << " ";
01323   }
01324       }
01325       os << "]" << endl;
01326 
01327       // mfh 16 Dec 2011: I found this commented out here; not sure if
01328       // we want to print this, so I'm leaving it commented out.
01329       /*
01330   os << "indices_from: ";
01331   k = 0;
01332   for( i = 0; i < nrecvs_; i++ )
01333   {
01334   for( j = 0; j < lengths_from_[i]; j++ )
01335   os << " " << indices_from_[j+k];
01336   k += lengths_from_[i];
01337   }
01338       */
01339 
01340       // Last output is a flush; it leaves a space and also 
01341       // helps synchronize output.
01342       os << std::flush;
01343     } // if it's my process' turn to print
01344 
01345     // Execute barriers to give output time to synchronize.
01346     // One barrier generally isn't enough.
01347     MPI_Barrier (comm_);
01348     MPI_Barrier (comm_);
01349     MPI_Barrier (comm_);
01350   }
01351 }
01352 
01353 //---------------------------------------------------------------------------
01354 int Epetra_MpiDistributor::Sort_ints_(
01355  int *vals_sort,     //  values to be sorted  
01356  int *vals_other,    // other array to be reordered with sort
01357  int  nvals)         // length of these two arrays
01358 {
01359 // It is primarily used to sort messages to improve communication flow. 
01360 // This routine will also insure that the ordering produced by the invert_map
01361 // routines is deterministic.  This should make bugs more reproducible.  This
01362 // is accomplished by sorting the message lists by processor ID.
01363 // This is a distribution count sort algorithm (see Knuth)
01364 //  This version assumes non negative integers. 
01365    
01366     if (nvals <= 1) return 0;
01367         
01368     int i;                        // loop counter        
01369      
01370     // find largest int, n, to size sorting array, then allocate and clear it
01371     int n = 0;  
01372     for (i = 0; i < nvals; i++)  
01373        if (n < vals_sort[i]) n = vals_sort[i]; 
01374     int *pos = new int [n+2];  
01375     for (i = 0; i < n+2; i++) pos[i] = 0;
01376  
01377     // copy input arrays into temporary copies to allow sorting original arrays
01378     int *copy_sort  = new int [nvals]; 
01379     int *copy_other = new int [nvals]; 
01380     for (i = 0; i < nvals; i++)
01381     { 
01382       copy_sort[i]  = vals_sort[i]; 
01383       copy_other[i] = vals_other[i];  
01384     }                           
01385  
01386     // count the occurances of integers ("distribution count")
01387     int *p = pos+1;
01388     for (i = 0; i < nvals; i++) p[copy_sort[i]]++;
01389  
01390     // create the partial sum of distribution counts 
01391     for (i = 1; i < n; i++) p[i] += p[i-1]; 
01392  
01393     // the shifted partitial sum is the index to store the data  in sort order
01394     p = pos; 
01395     for (i = 0; i < nvals; i++)         
01396     {                                   
01397       vals_sort  [p[copy_sort [i]]]   = copy_sort[i];
01398       vals_other [p[copy_sort [i]]++] = copy_other[i]; 
01399     } 
01400  
01401     delete [] copy_sort;
01402     delete [] copy_other; 
01403     delete [] pos; 
01404  
01405     return 0;
01406 }
01407 
01408 //-------------------------------------------------------------------------
01409 Epetra_MpiDistributor& Epetra_MpiDistributor::operator=(const Epetra_MpiDistributor& src)
01410 {
01411   (void)src;
01412   //not currently supported
01413   bool throw_error = true;
01414   if (throw_error) {
01415     throw ReportError("Epetra_MpiDistributor::operator= not supported.",-1);
01416   }
01417   return( *this );
01418 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines