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