Epetra_OffsetIndex.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 terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ************************************************************************
00027 //@HEADER
00028 
00029 #include "Epetra_OffsetIndex.h"
00030 #include "Epetra_CrsGraph.h"
00031 #include "Epetra_Import.h"
00032 #include "Epetra_Export.h"
00033 #include "Epetra_Distributor.h"
00034 #include "Epetra_Comm.h"
00035 
00036 //==============================================================================
00037 // Epetra_OffsetIndex constructor from Importer
00038 Epetra_OffsetIndex::Epetra_OffsetIndex( const Epetra_CrsGraph & SourceGraph,
00039                                         const Epetra_CrsGraph & TargetGraph,
00040                                         Epetra_Import & Importer )
00041   : Epetra_Object("Epetra::OffsetIndex"),
00042     NumSame_(0),
00043     SameOffsets_(0),
00044     NumPermute_(0),
00045     PermuteOffsets_(0),
00046     NumExport_(0),
00047     NumRemote_(0),
00048     RemoteOffsets_(0),
00049     DataOwned_(true)
00050 {
00051   NumSame_ = Importer.NumSameIDs();
00052 
00053   NumPermute_ = Importer.NumPermuteIDs();
00054   int * PermuteLIDs = Importer.PermuteToLIDs();
00055 
00056   NumExport_ = Importer.NumExportIDs();
00057   int * ExportLIDs = Importer.ExportLIDs();
00058 
00059   NumRemote_ = Importer.NumRemoteIDs();
00060   int * RemoteLIDs = Importer.RemoteLIDs();
00061 
00062   GenerateLocalOffsets_( SourceGraph, TargetGraph,
00063                          PermuteLIDs );
00064 
00065   GenerateRemoteOffsets_( SourceGraph, TargetGraph,
00066                           ExportLIDs, RemoteLIDs,
00067                           Importer.Distributor() );
00068 }
00069 
00070 //==============================================================================
00071 // Epetra_OffsetIndex constructor from Exporter
00072 Epetra_OffsetIndex::Epetra_OffsetIndex( const Epetra_CrsGraph & SourceGraph,
00073                                         const Epetra_CrsGraph & TargetGraph,
00074                                         Epetra_Export & Exporter )
00075   : Epetra_Object("Epetra::OffsetIndex"),
00076     NumSame_(0),
00077     SameOffsets_(0),
00078     NumPermute_(0),
00079     PermuteOffsets_(0),
00080     NumExport_(0),
00081     NumRemote_(0),
00082     RemoteOffsets_(0),
00083     DataOwned_(true)
00084 {
00085   NumSame_ = Exporter.NumSameIDs();
00086 
00087   NumPermute_ = Exporter.NumPermuteIDs();
00088   int * PermuteLIDs = Exporter.PermuteToLIDs();
00089 
00090   NumExport_ = Exporter.NumExportIDs();
00091   int * ExportLIDs = Exporter.ExportLIDs();
00092 
00093   NumRemote_ = Exporter.NumRemoteIDs();
00094   int * RemoteLIDs = Exporter.RemoteLIDs();
00095 
00096   GenerateLocalOffsets_( SourceGraph, TargetGraph,
00097                          PermuteLIDs );
00098 
00099   GenerateRemoteOffsets_( SourceGraph, TargetGraph,
00100                           ExportLIDs, RemoteLIDs,
00101                           Exporter.Distributor() );
00102 }
00103 
00104 //==============================================================================
00105 // Epetra_OffsetIndex copy constructor 
00106 Epetra_OffsetIndex::Epetra_OffsetIndex(const Epetra_OffsetIndex& Indexor)
00107   : Epetra_Object(Indexor),
00108     NumSame_(Indexor.NumSame_),
00109     SameOffsets_(Indexor.SameOffsets_),
00110     NumPermute_(Indexor.NumPermute_),
00111     PermuteOffsets_(Indexor.PermuteOffsets_),
00112     NumExport_(0),
00113     NumRemote_(Indexor.NumRemote_),
00114     RemoteOffsets_(Indexor.RemoteOffsets_),
00115     DataOwned_(false)
00116 {
00117 }
00118 
00119 //==============================================================================
00120 // Epetra_OffsetIndex destructor
00121 Epetra_OffsetIndex::~Epetra_OffsetIndex()
00122 {
00123   if( DataOwned_ )
00124   {
00125     for( int i = 0; i < NumSame_; ++i )
00126       if( SameOffsets_[i] ) delete [] SameOffsets_[i];
00127     delete [] SameOffsets_;
00128     for( int i = 0; i < NumPermute_; ++i )
00129       if( PermuteOffsets_[i] ) delete [] PermuteOffsets_[i];
00130     delete [] PermuteOffsets_;
00131     for( int i = 0; i < NumRemote_; ++i )
00132       if( RemoteOffsets_[i] ) delete [] RemoteOffsets_[i];
00133     delete [] RemoteOffsets_;
00134   }
00135 }
00136 
00137 //==============================================================================
00138 void Epetra_OffsetIndex::GenerateLocalOffsets_( const Epetra_CrsGraph & SourceGraph,
00139                                                 const Epetra_CrsGraph & TargetGraph,
00140                                                 const int * PermuteLIDs )
00141 {
00142   const int GlobalMaxNumSourceIndices = SourceGraph.GlobalMaxNumIndices();
00143 
00144   int NumSourceIndices;
00145   int * SourceIndices = 0;
00146   if( GlobalMaxNumSourceIndices>0 ) SourceIndices = new int[GlobalMaxNumSourceIndices];
00147 
00148   //setup Same Offsets
00149   SameOffsets_ = new int*[NumSame_];
00150   for( int i = 0; i < NumSame_; ++i ) SameOffsets_[i] = 0;
00151 
00152   for( int i = 0; i < NumSame_; ++i ) {
00153     int GID = SourceGraph.GRID(i);
00154     SourceGraph.ExtractGlobalRowCopy( GID,
00155                                       GlobalMaxNumSourceIndices,
00156                                       NumSourceIndices,
00157                                       SourceIndices );
00158 
00159     if( NumSourceIndices > 0 ) SameOffsets_[i] = new int[NumSourceIndices];
00160 
00161     int Loc = 0;
00162     int Start = 0;
00163     for( int j = 0; j < NumSourceIndices; ++j ) {
00164       Start = Loc;
00165       if( TargetGraph.FindGlobalIndexLoc(i,SourceIndices[j],Start,Loc) )
00166         SameOffsets_[i][j] = Loc;
00167       else
00168         SameOffsets_[i][j] = -1;
00169     }
00170   }
00171 
00172   //do also for permuted ids
00173   PermuteOffsets_ = new int*[NumPermute_];
00174   for( int i = 0; i < NumPermute_; ++i ) PermuteOffsets_[i] = 0;
00175 
00176   for( int i = 0; i < NumPermute_; ++i ) {
00177     int GID = SourceGraph.GRID(PermuteLIDs[i]);
00178     SourceGraph.ExtractGlobalRowCopy( GID,
00179                                       GlobalMaxNumSourceIndices,
00180                                       NumSourceIndices,
00181                                       SourceIndices );
00182 
00183     if( NumSourceIndices > 0 ) PermuteOffsets_[i] = new int[NumSourceIndices];
00184 
00185     int Loc = 0;
00186     int Start = 0;
00187     for( int j = 0; j < NumSourceIndices; ++j ) {
00188       Start = Loc;
00189       if( TargetGraph.FindGlobalIndexLoc(PermuteLIDs[i],SourceIndices[j],Start,Loc) )
00190         PermuteOffsets_[i][j] = Loc;
00191       else
00192         PermuteOffsets_[i][j] = -1;
00193     }
00194   }
00195 
00196   if( GlobalMaxNumSourceIndices>0 ) delete [] SourceIndices;
00197 }
00198 
00199 //==============================================================================
00200 void Epetra_OffsetIndex::GenerateRemoteOffsets_( const Epetra_CrsGraph & SourceGraph,
00201                                                  const Epetra_CrsGraph & TargetGraph,
00202                                                  const int * ExportLIDs,
00203                                                  const int * RemoteLIDs,
00204                                                  Epetra_Distributor & Distor )
00205 {
00206   int numProcs = SourceGraph.RowMap().Comm().NumProc();
00207   if (numProcs < 2) {
00208     return;
00209   }
00210 
00211   const int GlobalMaxNumIndices = SourceGraph.GlobalMaxNumIndices();
00212 
00213   int NumIndices;
00214   int * Indices = 0;
00215   if( GlobalMaxNumIndices>0 ) Indices = new int[GlobalMaxNumIndices];
00216 
00217   //Pack Source Rows
00218   int * Sizes = 0;
00219   if( NumExport_ > 0 ) Sizes = new int[NumExport_];
00220   int TotalSize = 0;
00221   for( int i = 0; i < NumExport_; ++i ) {
00222     Sizes[i] = SourceGraph.NumMyIndices(ExportLIDs[i]) + 1;
00223     TotalSize += Sizes[i];
00224   }
00225 
00226   int * SourceArray = new int[TotalSize+1];
00227   int Loc = 0;
00228   for( int i = 0; i < NumExport_; ++i ) {
00229     int GID = SourceGraph.GRID(ExportLIDs[i]);
00230     SourceArray[Loc] = Sizes[i]-1;
00231     SourceGraph.ExtractGlobalRowCopy( GID,
00232                                       GlobalMaxNumIndices,
00233                                       NumIndices,
00234                                       &(SourceArray[Loc+1]) );
00235     Loc += Sizes[i];
00236   }
00237 
00238   //Push to Target
00239   char * cRecvArray = 0;
00240   int * RecvArray = 0;
00241   int RecvArraySize = 0;
00242   Distor.Do( reinterpret_cast<char *>(SourceArray),
00243              (int)sizeof(int),
00244              Sizes,
00245              RecvArraySize,
00246              cRecvArray );
00247   RecvArray = reinterpret_cast<int*>(cRecvArray);
00248 
00249   //Construct RemoteOffsets
00250   if( NumRemote_ > 0 ) RemoteOffsets_ = new int*[NumRemote_];
00251   for( int i = 0; i < NumRemote_; ++i ) RemoteOffsets_[i] = 0;
00252 
00253   Loc = 0;
00254   for( int i = 0; i < NumRemote_; ++i ) {
00255     NumIndices = RecvArray[Loc];
00256     RemoteOffsets_[i] = new int[NumIndices];
00257     ++Loc;
00258     int FLoc = 0;
00259     int Start = 0;
00260     for( int j = 0; j < NumIndices; ++j ) {
00261       Start = FLoc;
00262       if( TargetGraph.FindGlobalIndexLoc(RemoteLIDs[i],RecvArray[Loc],Start,FLoc) )
00263         RemoteOffsets_[i][j] = FLoc;
00264       else
00265         RemoteOffsets_[i][j] = -1;
00266       ++Loc;
00267     }
00268   }
00269 
00270   if( GlobalMaxNumIndices>0 ) delete [] Indices;
00271   if( Sizes ) delete [] Sizes;
00272   if( SourceArray ) delete [] SourceArray;
00273   if( RecvArraySize ) delete [] cRecvArray;
00274 }
00275 
00276 //=============================================================================
00277 void Epetra_OffsetIndex::Print(ostream & os) const
00278 {
00279   os << "Number of Same IDs = " << NumSame_ << endl;
00280   os << "Number of Permute IDs = " << NumPermute_ << endl;
00281   os << "Number of Remote IDs = " << NumRemote_ << endl;
00282   
00283   return;
00284 }
00285 

Generated on Wed May 12 21:41:05 2010 for Epetra Package Browser (Single Doxygen Collection) by  doxygen 1.4.7