Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
test/MapColoring/cxx_main.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
43// Epetra_BlockMap Test routine
44
45#include "Epetra_BlockMap.h"
46#include "Epetra_Map.h"
47#include "Epetra_MapColoring.h"
48#include "Epetra_IntVector.h"
49#include "Epetra_Import.h"
50#ifdef EPETRA_MPI
51#include "Epetra_MpiComm.h"
52#include <mpi.h>
53#endif
54
55#include "Epetra_SerialComm.h"
56#include "Epetra_Time.h"
57#include "Epetra_Version.h"
58
59int main(int argc, char *argv[]) {
60
61 int i, returnierr=0;
62
63#ifdef EPETRA_MPI
64 // Initialize MPI
65 MPI_Init(&argc,&argv);
66 Epetra_MpiComm Comm(MPI_COMM_WORLD);
67#else
69#endif
70
71 // Uncomment to debug in parallel int tmp; if (Comm.MyPID()==0) cin >> tmp; Comm.Barrier();
72
73 bool verbose = false;
74 bool veryVerbose = false;
75
76 // Check if we should print results to standard out
77 if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
78
79 // Check if we should print lots of results to standard out
80 if (argc>2) if (argv[2][0]=='-' && argv[2][1]=='v') veryVerbose = true;
81
82 if (verbose && Comm.MyPID()==0)
83 std::cout << Epetra_Version() << std::endl << std::endl;
84
85 if (!verbose) Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
86
87 if (verbose) std::cout << Comm << std::endl << std::flush;
88
89 bool verbose1 = verbose;
90 if (verbose) verbose = (Comm.MyPID()==0);
91
92 bool veryVerbose1 = veryVerbose;
93 if (veryVerbose) veryVerbose = (Comm.MyPID()==0);
94
95 int NumMyElements = 100;
96 if (veryVerbose1) NumMyElements = 10;
97 NumMyElements += Comm.MyPID();
98 int MaxNumMyElements = NumMyElements+Comm.NumProc()-1;
99 int * ElementSizeList = new int[NumMyElements];
100 int * MyGlobalElements = new int[NumMyElements];
101
102 for (i = 0; i<NumMyElements; i++) {
103 MyGlobalElements[i] = (Comm.MyPID()*MaxNumMyElements+i)*2;
104 ElementSizeList[i] = i%6 + 2; // elementsizes go from 2 to 7
105 }
106
107 Epetra_BlockMap Map(-1, NumMyElements, MyGlobalElements, ElementSizeList,
108 0, Comm);
109
110 delete [] ElementSizeList;
111 delete [] MyGlobalElements;
112
113 Epetra_MapColoring C0(Map);
114
115 int * elementColors = new int[NumMyElements];
116
117 int maxcolor = 24;
118 int * colorCount = new int[maxcolor];
119 int ** colorLIDs = new int*[maxcolor];
120 for (i=0; i<maxcolor; i++) colorCount[i] = 0;
121 for (i=0; i<maxcolor; i++) colorLIDs[i] = 0;
122
123 for (i=0; i<Map.NumMyElements(); i++) {
124 assert(C0[i]==C0.DefaultColor());
125 assert(C0(Map.GID(i))==C0.DefaultColor());
126 if (i%2==0) C0[i] = i%6+5+i%14; // cycle through 5...23 on even elements
127 else C0(Map.GID(i)) = i%5+1; // cycle through 1...5 on odd elements
128 elementColors[i] = C0[i]; // Record color of ith element for use below
129 colorCount[C0[i]]++; // Count how many of each color for checking below
130 }
131
132 if (veryVerbose)
133 std::cout << "Original Map Coloring using element-by-element definitions" << std::endl;
134 if (veryVerbose1)
135 std::cout << C0 << std::endl;
136
137 int numColors = 0;
138 for (i=0; i<maxcolor; i++)
139 if (colorCount[i]>0) {
140 numColors++;
141 colorLIDs[i] = new int[colorCount[i]];
142 }
143 for (i=0; i<maxcolor; i++) colorCount[i] = 0;
144 for (i=0; i<Map.NumMyElements(); i++) colorLIDs[C0[i]][colorCount[C0[i]]++] = i;
145
146
147
148 int newDefaultColor = -1;
149 Epetra_MapColoring C1(Map, elementColors, newDefaultColor);
150 if (veryVerbose)
151 std::cout << "Same Map Coloring using one-time construction" << std::endl;
152 if (veryVerbose1)
153 std::cout << C1 << std::endl;
154 assert(C1.DefaultColor()==newDefaultColor);
155 for (i=0; i<Map.NumMyElements(); i++) assert(C1[i]==C0[i]);
156
157 Epetra_MapColoring C2(C1);
158 if (veryVerbose)
159 std::cout << "Same Map Coloring using copy constructor" << std::endl;
160 if (veryVerbose1)
161 std::cout << C1 << std::endl;
162 for (i=0; i<Map.NumMyElements(); i++) assert(C2[i]==C0[i]);
163 assert(C2.DefaultColor()==newDefaultColor);
164
165 assert(numColors==C2.NumColors());
166
167 for (i=0; i<maxcolor; i++) {
168 int curNumElementsWithColor = C2.NumElementsWithColor(i);
169 assert(colorCount[i]==curNumElementsWithColor);
170 if (curNumElementsWithColor==0) {
171 assert(C2.ColorLIDList(i)==0);
172 }
173 else
174 for (int j=0; j<curNumElementsWithColor; j++) assert(C2.ColorLIDList(i)[j]==colorLIDs[i][j]);
175 }
176 int curColor = 1;
177 Epetra_Map * Map1 = C2.GenerateMap(curColor);
178 Epetra_BlockMap * Map2 = C2.GenerateBlockMap(curColor);
179
180 assert(Map1->NumMyElements()==colorCount[curColor]);
181 assert(Map2->NumMyElements()==colorCount[curColor]);
182
183 for (i=0; i<Map1->NumMyElements(); i++) {
184 assert(Map1->GID(i)==Map.GID(colorLIDs[curColor][i]));
185 assert(Map2->GID(i)==Map.GID(colorLIDs[curColor][i]));
186 assert(Map2->ElementSize(i)==Map.ElementSize(colorLIDs[curColor][i]));
187 }
188
189 // Now test data redistribution capabilities
190
191
192 Epetra_Map ContiguousMap(-1, Map.NumMyElements(), Map.IndexBase(), Comm);
193 // This vector contains the element sizes for the original map.
194 Epetra_IntVector elementSizes(Copy, ContiguousMap, Map.ElementSizeList());
195 Epetra_IntVector elementIDs(Copy, ContiguousMap, Map.MyGlobalElements());
196 Epetra_IntVector elementColorValues(Copy, ContiguousMap, C2.ElementColors());
197
198
199 int NumMyElements0 = 0;
200 if (Comm.MyPID()==0) NumMyElements0 = Map.NumGlobalElements();
201 Epetra_Map CMap0(-1, NumMyElements0, Map.IndexBase(), Comm);
202 Epetra_Import importer(CMap0, ContiguousMap);
203 Epetra_IntVector elementSizes0(CMap0);
204 Epetra_IntVector elementIDs0(CMap0);
205 Epetra_IntVector elementColorValues0(CMap0);
206 elementSizes0.Import(elementSizes, importer, Insert);
207 elementIDs0.Import(elementIDs, importer, Insert);
208 elementColorValues0.Import(elementColorValues, importer, Insert);
209
210 Epetra_BlockMap MapOnPE0(-1,NumMyElements0, elementIDs0.Values(),
211 elementSizes0.Values(), Map.IndexBase(), Comm);
212
213 Epetra_Import importer1(MapOnPE0, Map);
214 Epetra_MapColoring ColoringOnPE0(MapOnPE0);
215 ColoringOnPE0.Import(C2, importer1, Insert);
216
217 for (i=0; i<MapOnPE0.NumMyElements(); i++)
218 assert(ColoringOnPE0[i]==elementColorValues0[i]);
219
220 if (veryVerbose)
221 std::cout << "Same Map Coloring on PE 0 only" << std::endl;
222 if (veryVerbose1)
223 std::cout << ColoringOnPE0 << std::endl;
224 Epetra_MapColoring C3(Map);
225 C3.Export(ColoringOnPE0, importer1, Insert);
226 for (i=0; i<Map.NumMyElements(); i++) assert(C3[i]==C2[i]);
227 if (veryVerbose)
228 std::cout << "Same Map Coloring after Import/Export exercise" << std::endl;
229 if (veryVerbose1)
230 std::cout << ColoringOnPE0 << std::endl;
231
232
233 if (verbose) std::cout << "Checked OK\n\n" << std::endl;
234
235 if (verbose1) {
236 if (verbose) std::cout << "Test ostream << operator" << std::endl << std::flush;
237 std::cout << C0 << std::endl;
238 }
239
240
241 delete [] elementColors;
242 for (i=0; i<maxcolor; i++) if (colorLIDs[i]!=0) delete [] colorLIDs[i];
243 delete [] colorLIDs;
244 delete [] colorCount;
245
246 delete Map1;
247 delete Map2;
248
249
250#ifdef EPETRA_MPI
251 MPI_Finalize();
252#endif
253
254 return returnierr;
255}
256
std::string Epetra_Version()
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
int * ElementSizeList() const
List of the element sizes corresponding to the array MyGlobalElements().
int GID(int LID) const
Returns global ID of local ID, return IndexBase-1 if not found on this processor.
int IndexBase() const
Index base for this map.
int ElementSize() const
Returns the size of elements in the map; only valid if map has constant element size.
int NumGlobalElements() const
Number of elements across all processors.
int NumMyElements() const
Number of elements on the calling processor.
int Import(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Imports an Epetra_DistObject using the Epetra_Import object.
int Export(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Exports an Epetra_DistObject using the Epetra_Import object.
Epetra_Import: This class builds an import object for efficient importing of off-processor elements.
Epetra_IntVector: A class for constructing and using dense integer vectors on a parallel computer.
int * Values() const
Returns a pointer to an array containing the values of this vector.
Epetra_MapColoring: A class for coloring Epetra_Map and Epetra_BlockMap objects.
int NumColors() const
Returns number of colors on the calling processor.
int NumElementsWithColor(int Color) const
Returns number of map elements on calling processor having specified Color.
int DefaultColor() const
Returns default color.
Epetra_BlockMap * GenerateBlockMap(int Color) const
Generates an Epetra_BlockMap of the GIDs associated with the specified color.
Epetra_Map * GenerateMap(int Color) const
Generates an Epetra_Map of the GIDs associated with the specified color.
int * ElementColors() const
Returns pointer to array of the colors associated with the LIDs on the calling processor.
int * ColorLIDList(int Color) const
Returns pointer to array of Map LIDs associated with the specified color.
Epetra_Map: A class for partitioning vectors and matrices.
Definition Epetra_Map.h:119
Epetra_MpiComm: The Epetra MPI Communication Class.
int NumProc() const
Returns total number of processes.
int MyPID() const
Return my process ID.
static void SetTracebackMode(int TracebackModeValue)
Set the value of the Epetra_Object error traceback report mode.
Epetra_SerialComm: The Epetra Serial Communication Class.
int main(int argc, char *argv[])