test/MapColoring/cxx_main.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 // Epetra_BlockMap Test routine
00030 
00031 #include "Epetra_BlockMap.h"
00032 #include "Epetra_Map.h"
00033 #include "Epetra_MapColoring.h"
00034 #include "Epetra_IntVector.h"
00035 #include "Epetra_Import.h"
00036 #ifdef EPETRA_MPI
00037 #include "Epetra_MpiComm.h"
00038 #include <mpi.h>
00039 #endif
00040 
00041 #include "Epetra_SerialComm.h"
00042 #include "Epetra_Time.h"
00043 #include "Epetra_Version.h"
00044 
00045 int main(int argc, char *argv[]) {
00046 
00047   int i, returnierr=0;
00048 
00049 #ifdef EPETRA_MPI
00050   // Initialize MPI
00051   MPI_Init(&argc,&argv);
00052   Epetra_MpiComm Comm(MPI_COMM_WORLD);
00053 #else
00054   Epetra_SerialComm Comm;
00055 #endif
00056 
00057   // Uncomment to debug in parallel int tmp; if (Comm.MyPID()==0) cin >> tmp; Comm.Barrier();
00058 
00059   bool verbose = false;
00060   bool veryVerbose = false;
00061 
00062   // Check if we should print results to standard out
00063   if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
00064 
00065   // Check if we should print lots of results to standard out
00066   if (argc>2) if (argv[2][0]=='-' && argv[2][1]=='v') veryVerbose = true;
00067 
00068   if (verbose && Comm.MyPID()==0)
00069     cout << Epetra_Version() << endl << endl;
00070 
00071   if (!verbose) Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
00072 
00073   if (verbose) cout << Comm << endl << flush;
00074 
00075   bool verbose1 = verbose;
00076   if (verbose) verbose = (Comm.MyPID()==0);
00077 
00078   bool veryVerbose1 = veryVerbose;
00079   if (veryVerbose) veryVerbose = (Comm.MyPID()==0);
00080 
00081   int NumMyElements = 100;
00082   if (veryVerbose1) NumMyElements = 10;
00083   NumMyElements += Comm.MyPID();
00084   int MaxNumMyElements = NumMyElements+Comm.NumProc()-1;
00085   int * ElementSizeList = new int[NumMyElements];
00086   int * MyGlobalElements = new int[NumMyElements];
00087 
00088   for (i = 0; i<NumMyElements; i++) {
00089     MyGlobalElements[i] = (Comm.MyPID()*MaxNumMyElements+i)*2;
00090     ElementSizeList[i] = i%6 + 2; // elementsizes go from 2 to 7
00091   }
00092 
00093   Epetra_BlockMap Map(-1, NumMyElements, MyGlobalElements, ElementSizeList,
00094           0, Comm);
00095 
00096   delete [] ElementSizeList;
00097   delete [] MyGlobalElements;
00098 
00099   Epetra_MapColoring C0(Map);
00100 
00101   int * elementColors = new int[NumMyElements];
00102 
00103   int maxcolor = 24;
00104   int * colorCount = new int[maxcolor];
00105   int ** colorLIDs = new int*[maxcolor];
00106   for (i=0; i<maxcolor; i++) colorCount[i] = 0;
00107   for (i=0; i<maxcolor; i++) colorLIDs[i] = 0;
00108 
00109   int defaultColor = C0.DefaultColor();
00110   for (i=0; i<Map.NumMyElements(); i++) {
00111     assert(C0[i]==defaultColor);
00112     assert(C0(Map.GID(i))==defaultColor);
00113     if (i%2==0) C0[i] = i%6+5+i%14; // cycle through 5...23 on even elements
00114     else C0(Map.GID(i)) = i%5+1; // cycle through 1...5 on odd elements
00115     elementColors[i] = C0[i]; // Record color of ith element for use below
00116     colorCount[C0[i]]++; // Count how many of each color for checking below
00117   }
00118   
00119   if (veryVerbose)
00120     cout << "Original Map Coloring using element-by-element definitions" << endl;
00121   if (veryVerbose1)
00122     cout <<  C0 << endl;
00123 
00124   int numColors = 0;
00125   for (i=0; i<maxcolor; i++) 
00126     if (colorCount[i]>0) {
00127       numColors++;
00128       colorLIDs[i] = new int[colorCount[i]];
00129     }
00130   for (i=0; i<maxcolor; i++) colorCount[i] = 0;
00131   for (i=0; i<Map.NumMyElements(); i++) colorLIDs[C0[i]][colorCount[C0[i]]++] = i;
00132 
00133   
00134 
00135   int newDefaultColor = -1;
00136   Epetra_MapColoring C1(Map, elementColors, newDefaultColor);
00137   if (veryVerbose)
00138     cout << "Same Map Coloring using one-time construction" << endl;
00139   if (veryVerbose1)
00140     cout <<  C1 << endl;
00141   assert(C1.DefaultColor()==newDefaultColor);
00142   for (i=0; i<Map.NumMyElements(); i++) assert(C1[i]==C0[i]);
00143 
00144   Epetra_MapColoring C2(C1);
00145   if (veryVerbose)
00146     cout << "Same Map Coloring using copy constructor" << endl;
00147   if (veryVerbose1)
00148     cout <<  C1 << endl;
00149   for (i=0; i<Map.NumMyElements(); i++) assert(C2[i]==C0[i]);
00150   assert(C2.DefaultColor()==newDefaultColor);
00151 
00152   assert(numColors==C2.NumColors());
00153 
00154   for (i=0; i<maxcolor; i++) {
00155     int curNumElementsWithColor = C2.NumElementsWithColor(i);
00156     assert(colorCount[i]==curNumElementsWithColor);
00157     int * curColorLIDList = C2.ColorLIDList(i);
00158     if (curNumElementsWithColor==0) {
00159       assert(curColorLIDList==0);
00160     }
00161     else
00162       for (int j=0; j<curNumElementsWithColor; j++) assert(curColorLIDList[j]==colorLIDs[i][j]);
00163   }
00164   int curColor = 1;
00165   Epetra_Map * Map1 = C2.GenerateMap(curColor);
00166   Epetra_BlockMap * Map2 = C2.GenerateBlockMap(curColor);
00167 
00168   assert(Map1->NumMyElements()==colorCount[curColor]);
00169   assert(Map2->NumMyElements()==colorCount[curColor]);
00170 
00171   for (i=0; i<Map1->NumMyElements(); i++) {
00172     assert(Map1->GID(i)==Map.GID(colorLIDs[curColor][i]));
00173     assert(Map2->GID(i)==Map.GID(colorLIDs[curColor][i]));
00174     assert(Map2->ElementSize(i)==Map.ElementSize(colorLIDs[curColor][i]));
00175   }
00176 
00177   // Now test data redistribution capabilities
00178 
00179 
00180   Epetra_Map ContiguousMap(-1, Map.NumMyElements(), Map.IndexBase(), Comm);
00181   // This vector contains the element sizes for the original map.
00182   Epetra_IntVector elementSizes(Copy, ContiguousMap, Map.ElementSizeList());
00183   Epetra_IntVector elementIDs(Copy, ContiguousMap, Map.MyGlobalElements());
00184   Epetra_IntVector elementColorValues(Copy, ContiguousMap, C2.ElementColors());
00185 
00186 
00187   int NumMyElements0 = 0;
00188   if (Comm.MyPID()==0) NumMyElements0 = Map.NumGlobalElements();
00189   Epetra_Map CMap0(-1, NumMyElements0, Map.IndexBase(), Comm);
00190   Epetra_Import importer(CMap0, ContiguousMap);
00191   Epetra_IntVector elementSizes0(CMap0);
00192   Epetra_IntVector elementIDs0(CMap0);
00193   Epetra_IntVector elementColorValues0(CMap0);
00194   elementSizes0.Import(elementSizes, importer, Insert);
00195   elementIDs0.Import(elementIDs, importer, Insert);
00196   elementColorValues0.Import(elementColorValues, importer, Insert);
00197 
00198   Epetra_BlockMap MapOnPE0(-1,NumMyElements0, elementIDs0.Values(), 
00199          elementSizes0.Values(), Map.IndexBase(), Comm);
00200 
00201   Epetra_Import importer1(MapOnPE0, Map);
00202   Epetra_MapColoring ColoringOnPE0(MapOnPE0);
00203   ColoringOnPE0.Import(C2, importer1, Insert);
00204 
00205   for (i=0; i<MapOnPE0.NumMyElements(); i++)
00206     assert(ColoringOnPE0[i]==elementColorValues0[i]);
00207 
00208   if (veryVerbose)
00209     cout << "Same Map Coloring on PE 0 only" << endl;
00210   if (veryVerbose1)
00211     cout <<  ColoringOnPE0 << endl;
00212   Epetra_MapColoring C3(Map);
00213   C3.Export(ColoringOnPE0, importer1, Insert);
00214   for (i=0; i<Map.NumMyElements(); i++) assert(C3[i]==C2[i]);
00215   if (veryVerbose)
00216     cout << "Same Map Coloring after Import/Export exercise" << endl;
00217   if (veryVerbose1)
00218     cout <<  ColoringOnPE0 << endl;
00219    
00220   
00221   if (verbose) cout << "Checked OK\n\n" <<endl;
00222 
00223   if (verbose1) {
00224     if (verbose) cout << "Test ostream << operator" << endl << flush;
00225     cout << C0 << endl;
00226   }
00227   
00228 
00229   delete [] elementColors;
00230   for (i=0; i<maxcolor; i++) if (colorLIDs[i]!=0) delete [] colorLIDs[i];
00231   delete [] colorLIDs;
00232   delete [] colorCount;
00233 
00234   delete Map1;
00235   delete Map2;
00236 
00237 
00238 #ifdef EPETRA_MPI
00239   MPI_Finalize();
00240 #endif
00241 
00242   return returnierr;
00243 }
00244 

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