Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Epetra_OffsetIndex.cpp
Go to the documentation of this file.
1//@HEADER
2// ************************************************************************
3//
4// Epetra: Linear Algebra Services Package
5// Copyright 2011 Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ************************************************************************
40//@HEADER
41
42#include "Epetra_ConfigDefs.h"
43#include "Epetra_OffsetIndex.h"
44#include "Epetra_CrsGraph.h"
45#include "Epetra_Import.h"
46#include "Epetra_Export.h"
47#include "Epetra_Distributor.h"
48#include "Epetra_Comm.h"
49
50//==============================================================================
51// Epetra_OffsetIndex constructor from Importer
53 const Epetra_CrsGraph & TargetGraph,
54 Epetra_Import & Importer )
55 : Epetra_Object("Epetra::OffsetIndex"),
56 NumSame_(0),
57 SameOffsets_(0),
58 NumPermute_(0),
59 PermuteOffsets_(0),
60 NumExport_(0),
61 NumRemote_(0),
62 RemoteOffsets_(0),
63 DataOwned_(true)
64{
65 NumSame_ = Importer.NumSameIDs();
66
67 NumPermute_ = Importer.NumPermuteIDs();
68 int * PermuteLIDs = Importer.PermuteToLIDs();
69
70 NumExport_ = Importer.NumExportIDs();
71 int * ExportLIDs = Importer.ExportLIDs();
72
73 NumRemote_ = Importer.NumRemoteIDs();
74 int * RemoteLIDs = Importer.RemoteLIDs();
75
76 if(!SourceGraph.RowMap().GlobalIndicesTypeMatch(TargetGraph.RowMap()))
77 throw ReportError("Epetra_OffsetIndex::Epetra_OffsetIndex: SourceGraph and TargetGraph global indices type mismatch", -1);
78 if(SourceGraph.RowMap().GlobalIndicesInt()) {
79#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
80 GenerateLocalOffsets_<int>( SourceGraph, TargetGraph,
81 PermuteLIDs );
82
83 GenerateRemoteOffsets_<int>( SourceGraph, TargetGraph,
84 ExportLIDs, RemoteLIDs,
85 Importer.Distributor() );
86#else
87 throw ReportError("Epetra_OffsetIndex::Epetra_OffsetIndex: ERROR, GlobalIndicesInt but no API for it.",-1);
88#endif
89 }
90 else if(SourceGraph.RowMap().GlobalIndicesLongLong()) {
91#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
92 GenerateLocalOffsets_<long long>( SourceGraph, TargetGraph,
93 PermuteLIDs );
94
95 GenerateRemoteOffsets_<long long>( SourceGraph, TargetGraph,
96 ExportLIDs, RemoteLIDs,
97 Importer.Distributor() );
98#else
99 throw ReportError("Epetra_OffsetIndex::Epetra_OffsetIndex: ERROR, GlobalIndicesLongLong but no API for it.",-1);
100#endif
101 }
102 else
103 throw ReportError("Epetra_OffsetIndex::Epetra_OffsetIndex: SourceGraph global indices type unknown", -1);
104}
105
106//==============================================================================
107// Epetra_OffsetIndex constructor from Exporter
109 const Epetra_CrsGraph & TargetGraph,
110 Epetra_Export & Exporter )
111 : Epetra_Object("Epetra::OffsetIndex"),
112 NumSame_(0),
113 SameOffsets_(0),
114 NumPermute_(0),
115 PermuteOffsets_(0),
116 NumExport_(0),
117 NumRemote_(0),
118 RemoteOffsets_(0),
119 DataOwned_(true)
120{
121 NumSame_ = Exporter.NumSameIDs();
122
123 NumPermute_ = Exporter.NumPermuteIDs();
124 int * PermuteLIDs = Exporter.PermuteToLIDs();
125
126 NumExport_ = Exporter.NumExportIDs();
127 int * ExportLIDs = Exporter.ExportLIDs();
128
129 NumRemote_ = Exporter.NumRemoteIDs();
130 int * RemoteLIDs = Exporter.RemoteLIDs();
131
132 if(!SourceGraph.RowMap().GlobalIndicesTypeMatch(TargetGraph.RowMap()))
133 throw ReportError("Epetra_OffsetIndex::Epetra_OffsetIndex: SourceGraph and TargetGraph global indices type mismatch", -1);
134 if(SourceGraph.RowMap().GlobalIndicesInt()) {
135#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
136 GenerateLocalOffsets_<int>( SourceGraph, TargetGraph,
137 PermuteLIDs );
138
139 GenerateRemoteOffsets_<int>( SourceGraph, TargetGraph,
140 ExportLIDs, RemoteLIDs,
141 Exporter.Distributor() );
142#else
143 throw ReportError("Epetra_OffsetIndex::Epetra_OffsetIndex: ERROR, GlobalIndicesInt but no API for it.",-1);
144#endif
145 }
146 else if(SourceGraph.RowMap().GlobalIndicesLongLong()) {
147#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
148 GenerateLocalOffsets_<long long>( SourceGraph, TargetGraph,
149 PermuteLIDs );
150
151 GenerateRemoteOffsets_<long long>( SourceGraph, TargetGraph,
152 ExportLIDs, RemoteLIDs,
153 Exporter.Distributor() );
154#else
155 throw ReportError("Epetra_OffsetIndex::Epetra_OffsetIndex: ERROR, GlobalIndicesLongLong but no API for it.",-1);
156#endif
157 }
158 else
159 throw ReportError("Epetra_OffsetIndex::Epetra_OffsetIndex: SourceGraph global indices type unknown", -1);
160}
161
162//==============================================================================
163// Epetra_OffsetIndex copy constructor
165 : Epetra_Object(Indexor),
166 NumSame_(Indexor.NumSame_),
167 SameOffsets_(Indexor.SameOffsets_),
168 NumPermute_(Indexor.NumPermute_),
169 PermuteOffsets_(Indexor.PermuteOffsets_),
170 NumExport_(0),
171 NumRemote_(Indexor.NumRemote_),
172 RemoteOffsets_(Indexor.RemoteOffsets_),
173 DataOwned_(false)
174{
175}
176
177//==============================================================================
178// Epetra_OffsetIndex destructor
180{
181 if( DataOwned_ )
182 {
183 for( int i = 0; i < NumSame_; ++i )
184 if( SameOffsets_[i] ) delete [] SameOffsets_[i];
185 delete [] SameOffsets_;
186 for( int i = 0; i < NumPermute_; ++i )
187 if( PermuteOffsets_[i] ) delete [] PermuteOffsets_[i];
188 delete [] PermuteOffsets_;
189 for( int i = 0; i < NumRemote_; ++i )
190 if( RemoteOffsets_[i] ) delete [] RemoteOffsets_[i];
191 delete [] RemoteOffsets_;
192 }
193}
194
195//==============================================================================
196template<typename int_type>
198 const Epetra_CrsGraph & TargetGraph,
199 const int * PermuteLIDs )
200{
201 const int GlobalMaxNumSourceIndices = SourceGraph.GlobalMaxNumIndices();
202
203 int NumSourceIndices;
204 int_type * SourceIndices = 0;
205 if( GlobalMaxNumSourceIndices>0 ) SourceIndices = new int_type[GlobalMaxNumSourceIndices];
206
207 //setup Same Offsets
208 SameOffsets_ = new int*[NumSame_];
209 for( int i = 0; i < NumSame_; ++i ) SameOffsets_[i] = 0;
210
211 for( int i = 0; i < NumSame_; ++i ) {
212 int_type GID = (int_type) SourceGraph.GRID64(i);
213 SourceGraph.ExtractGlobalRowCopy( GID,
214 GlobalMaxNumSourceIndices,
215 NumSourceIndices,
216 SourceIndices );
217
218 if( NumSourceIndices > 0 ) SameOffsets_[i] = new int[NumSourceIndices];
219
220 int Loc = 0;
221 int Start = 0;
222 for( int j = 0; j < NumSourceIndices; ++j ) {
223 Start = Loc;
224 if( TargetGraph.FindGlobalIndexLoc(i,SourceIndices[j],Start,Loc) )
225 SameOffsets_[i][j] = Loc;
226 else
227 SameOffsets_[i][j] = -1;
228 }
229 }
230
231 //do also for permuted ids
232 PermuteOffsets_ = new int*[NumPermute_];
233 for( int i = 0; i < NumPermute_; ++i ) PermuteOffsets_[i] = 0;
234
235 for( int i = 0; i < NumPermute_; ++i ) {
236 int_type GID = (int_type) SourceGraph.GRID64(PermuteLIDs[i]);
237 SourceGraph.ExtractGlobalRowCopy( GID,
238 GlobalMaxNumSourceIndices,
239 NumSourceIndices,
240 SourceIndices );
241
242 if( NumSourceIndices > 0 ) PermuteOffsets_[i] = new int[NumSourceIndices];
243
244 int Loc = 0;
245 int Start = 0;
246 for( int j = 0; j < NumSourceIndices; ++j ) {
247 Start = Loc;
248 if( TargetGraph.FindGlobalIndexLoc(PermuteLIDs[i],SourceIndices[j],Start,Loc) )
249 PermuteOffsets_[i][j] = Loc;
250 else
251 PermuteOffsets_[i][j] = -1;
252 }
253 }
254
255 if( GlobalMaxNumSourceIndices>0 ) delete [] SourceIndices;
256}
257
258//==============================================================================
259template<typename int_type>
261 const Epetra_CrsGraph & TargetGraph,
262 const int * ExportLIDs,
263 const int * RemoteLIDs,
264 Epetra_Distributor & Distor )
265{
266 int numProcs = SourceGraph.RowMap().Comm().NumProc();
267 if (numProcs < 2) {
268 return;
269 }
270
271 const int GlobalMaxNumIndices = SourceGraph.GlobalMaxNumIndices();
272
273 int NumIndices;
274 /* "Indices" appears to be unused -- jhurani@txcorp.com
275 int * Indices = 0;
276 if( GlobalMaxNumIndices>0 ) Indices = new int[GlobalMaxNumIndices];
277 */
278
279 //Pack Source Rows
280 int * Sizes = 0;
281 if( NumExport_ > 0 ) Sizes = new int[NumExport_];
282 int TotalSize = 0;
283 for( int i = 0; i < NumExport_; ++i ) {
284 Sizes[i] = SourceGraph.NumMyIndices(ExportLIDs[i]) + 1;
285 TotalSize += Sizes[i];
286 }
287
288 int_type * SourceArray = new int_type[TotalSize+1];
289 int Loc = 0;
290 for( int i = 0; i < NumExport_; ++i ) {
291 int_type GID = (int_type) SourceGraph.GRID64(ExportLIDs[i]);
292 SourceArray[Loc] = Sizes[i]-1;
293 SourceGraph.ExtractGlobalRowCopy( GID,
294 GlobalMaxNumIndices,
295 NumIndices,
296 &(SourceArray[Loc+1]) );
297 Loc += Sizes[i];
298 }
299
300 //Push to Target
301 char * cRecvArray = 0;
302 int_type * RecvArray = 0;
303 int RecvArraySize = 0;
304
305 Distor.Do( reinterpret_cast<char *>(SourceArray),
306 (int)sizeof(int_type),
307 Sizes,
308 RecvArraySize,
309 cRecvArray );
310 RecvArray = reinterpret_cast<int_type*>(cRecvArray);
311
312 //Construct RemoteOffsets
313 if( NumRemote_ > 0 ) RemoteOffsets_ = new int*[NumRemote_];
314 for( int i = 0; i < NumRemote_; ++i ) RemoteOffsets_[i] = 0;
315
316 Loc = 0;
317 for( int i = 0; i < NumRemote_; ++i ) {
318 NumIndices = (int) RecvArray[Loc];
319 RemoteOffsets_[i] = new int[NumIndices];
320 ++Loc;
321 int FLoc = 0;
322 int Start = 0;
323 for( int j = 0; j < NumIndices; ++j ) {
324 Start = FLoc;
325 if( TargetGraph.FindGlobalIndexLoc(RemoteLIDs[i],RecvArray[Loc],Start,FLoc) )
326 RemoteOffsets_[i][j] = FLoc;
327 else
328 RemoteOffsets_[i][j] = -1;
329 ++Loc;
330 }
331 }
332
333 /* "Indices" appears to be unused -- jhurani@txcorp.com
334 if( GlobalMaxNumIndices>0 ) delete [] Indices;
335 */
336
337 if( Sizes ) delete [] Sizes;
338 if( SourceArray ) delete [] SourceArray;
339 if( RecvArraySize ) delete [] cRecvArray;
340}
341
342//=============================================================================
343void Epetra_OffsetIndex::Print(std::ostream & os) const
344{
345 os << "Number of Same IDs = " << NumSame_ << std::endl;
346 os << "Number of Permute IDs = " << NumPermute_ << std::endl;
347 os << "Number of Remote IDs = " << NumRemote_ << std::endl;
348
349 return;
350}
351
bool GlobalIndicesTypeMatch(const Epetra_BlockMap &other) const
bool GlobalIndicesInt() const
Returns true if map create with int NumGlobalElements.
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
bool GlobalIndicesLongLong() const
Returns true if map create with long long NumGlobalElements.
virtual int NumProc() const =0
Returns total number of processes.
Epetra_CrsGraph: A class for constructing and using sparse compressed row graphs.
const Epetra_BlockMap & RowMap() const
Returns the RowMap associated with this graph.
bool FindGlobalIndexLoc(int LocalRow, int Index, int Start, int &Loc) const
int GlobalMaxNumIndices() const
Returns the maximun number of nonzero entries across all rows across all processors.
int NumMyIndices(int Row) const
Returns the current number of nonzero entries in specified local row on this processor.
int ExtractGlobalRowCopy(int GlobalRow, int LenOfIndices, int &NumIndices, int *Indices) const
Extract a list of elements in a specified global row of the graph. Put into storage allocated by call...
long long GRID64(int LRID_in) const
Epetra_Distributor: The Epetra Gather/Scatter Setup Base Class.
virtual int Do(char *export_objs, int obj_size, int &len_import_objs, char *&import_objs)=0
Execute plan on buffer of export objects in a single step.
Epetra_Export: This class builds an export object for efficient exporting of off-processor elements.
Definition: Epetra_Export.h:62
int * RemoteLIDs() const
List of elements in the target map that are coming from other processors.
Epetra_Distributor & Distributor() const
int NumRemoteIDs() const
Returns the number of elements that are not on the calling processor.
int * PermuteToLIDs() const
List of elements in the target map that are permuted.
int NumPermuteIDs() const
Returns the number of elements that are local to the calling processor, but not part of the first Num...
int NumSameIDs() const
Returns the number of elements that are identical between the source and target maps,...
int * ExportLIDs() const
List of elements that will be sent to other processors.
int NumExportIDs() const
Returns the number of elements that must be sent by the calling processor to other processors.
Epetra_Import: This class builds an import object for efficient importing of off-processor elements.
Definition: Epetra_Import.h:63
int * PermuteToLIDs() const
List of elements in the target map that are permuted.
int NumPermuteIDs() const
Returns the number of elements that are local to the calling processor, but not part of the first Num...
int NumExportIDs() const
Returns the number of elements that must be sent by the calling processor to other processors.
int NumSameIDs() const
Returns the number of elements that are identical between the source and target maps,...
Epetra_Distributor & Distributor() const
int * RemoteLIDs() const
List of elements in the target map that are coming from other processors.
int NumRemoteIDs() const
Returns the number of elements that are not on the calling processor.
int * ExportLIDs() const
List of elements that will be sent to other processors.
Epetra_Object: The base Epetra class.
Definition: Epetra_Object.h:57
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
Epetra_OffsetIndex: This class builds index for efficient mapping of data from one Epetra_CrsGraph ba...
void GenerateLocalOffsets_(const Epetra_CrsGraph &SourceGraph, const Epetra_CrsGraph &TargetGraph, const int *PermuteLIDs)
void GenerateRemoteOffsets_(const Epetra_CrsGraph &SourceGraph, const Epetra_CrsGraph &TargetGraph, const int *ExportLIDs, const int *RemoteLIDs, Epetra_Distributor &Distor)
Epetra_OffsetIndex(const Epetra_CrsGraph &SourceGraph, const Epetra_CrsGraph &TargetGraph, Epetra_Import &Importer)
Constructs a Epetra_OffsetIndex object from the graphs and an importer.
virtual void Print(std::ostream &os) const
Print object to an output stream Print method.
virtual ~Epetra_OffsetIndex(void)
Epetra_OffsetIndex destructor.