Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
rcbTest.cpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Zoltan2: A package of combinatorial algorithms for scientific computing
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Karen Devine (kddevin@sandia.gov)
39 // Erik Boman (egboman@sandia.gov)
40 // Siva Rajamanickam (srajama@sandia.gov)
41 //
42 // ***********************************************************************
43 //
44 // @HEADER
45 
51 #include <Zoltan2_config.h>
52 #include <Zoltan2_TestHelpers.hpp>
57 
58 using Teuchos::RCP;
59 using Teuchos::rcp;
60 
61 typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> tMVector_t;
63 // Need to use zgno_t as zzgid_t in basic types, since zzgid_t=zgno_t in MultiVector
64 
65 
74  const RCP<const Teuchos::Comm<int> > & comm,
75  int nParts,
76  string &filename,
77  bool doRemap
78 )
79 {
80  int me = comm->getRank();
81  if (me == 0)
82  std::cout << "Parallel partitioning of " << filename << ".mtx: "
83  << nParts << " parts." << std::endl;
84 
85  std::string fname(filename);
86  UserInputForTests uinput(testDataFilePath, fname, comm, true);
87 
88  RCP<tMVector_t> coords = uinput.getUICoordinates();
89  if (me == 0)
90  std::cout << "Multivector length = " << coords->getGlobalLength()
91  << " Num vectors = " << coords->getNumVectors() << std::endl;
92 
93  RCP<const tMVector_t> coordsConst = rcp_const_cast<const tMVector_t>(coords);
94 
95  typedef Zoltan2::XpetraMultiVectorAdapter<tMVector_t> inputAdapter_t;
96  inputAdapter_t ia(coordsConst);
97  if (me == 0)
98  std::cout << "Adapter constructed" << std::endl;
99 
100  Teuchos::ParameterList params("test params");
101  params.set("debug_level", "basic_status");
102  params.set("num_global_parts", nParts);
103  params.set("algorithm", "rcb");
104  params.set("imbalance_tolerance", 1.1);
105  if (doRemap) params.set("remap_parts", true); // bool parameter
106 
107 #ifdef HAVE_ZOLTAN2_MPI
109  MPI_COMM_WORLD);
110 #else
111  Zoltan2::PartitioningProblem<inputAdapter_t> problem(&ia, &params);
112 #endif
113  if (me == 0)
114  std::cout << "Problem constructed" << std::endl;
115 
116 
117  problem.solve();
118  if (me == 0)
119  std::cout << "Problem solved" << std::endl;
120 }
121 
122 void serialTest(int numParts, bool doRemap)
123 {
124  int numCoords = 1000;
125  numParts *= 8;
126 
127  std::cout << "Serial partitioning: " << numParts << " parts." << std::endl;
128 
129  zgno_t *ids = new zgno_t [numCoords];
130  if (!ids)
131  throw std::bad_alloc();
132  for (int i=0; i < numCoords; i++)
133  ids[i] = i;
134  ArrayRCP<zgno_t> globalIds(ids, 0, numCoords, true);
135 
136  Array<ArrayRCP<zscalar_t> > randomCoords(3);
137  UserInputForTests::getUIRandomData(555, numCoords, 0, 10,
138  randomCoords.view(0,3));
139 
140  typedef Zoltan2::BasicVectorAdapter<myTypes_t> inputAdapter_t;
141 
142  inputAdapter_t ia(numCoords, ids,
143  randomCoords[0].getRawPtr(), randomCoords[1].getRawPtr(),
144  randomCoords[2].getRawPtr(), 1,1,1);
145 
146  Teuchos::ParameterList params("test params");
147  params.set("debug_level", "basic_status");
148  params.set("num_global_parts", numParts);
149  params.set("algorithm", "rcb");
150  params.set("imbalance_tolerance", 1.1);
151  if (doRemap) params.set("remap_parts", true); // bool parameter
152 
153 #ifdef HAVE_ZOLTAN2_MPI
155  &ia, &params, MPI_COMM_SELF);
156 #else
157  Zoltan2::PartitioningProblem<inputAdapter_t> serialProblem(&ia, &params);
158 #endif
159 
160  serialProblem.solve();
161 }
162 
163 void meshCoordinatesTest(const RCP<const Teuchos::Comm<int> > & comm)
164 {
165  int xdim = 40;
166  int ydim = 60;
167  int zdim = 20;
168  UserInputForTests uinput(xdim, ydim, zdim, string("Laplace3D"), comm, true, true);
169 
170  RCP<tMVector_t> coords = uinput.getUICoordinates();
171 
172  size_t localCount = coords->getLocalLength();
173 
174  zscalar_t *x=NULL, *y=NULL, *z=NULL;
175  x = coords->getDataNonConst(0).getRawPtr();
176  y = coords->getDataNonConst(1).getRawPtr();
177  z = coords->getDataNonConst(2).getRawPtr();
178 
179  const zgno_t *globalIds = coords->getMap()->getNodeElementList().getRawPtr();
180  typedef Zoltan2::BasicVectorAdapter<tMVector_t> inputAdapter_t;
181 
182  inputAdapter_t ia(localCount, globalIds, x, y, z, 1, 1, 1);
183 
184  Teuchos::ParameterList params("test params");
185  params.set("rectilinear", true); // bool parameter
186 
187 #ifdef HAVE_ZOLTAN2_MPI
188  Zoltan2::PartitioningProblem<inputAdapter_t> problem(&ia, &params, MPI_COMM_WORLD);
189 #else
190  Zoltan2::PartitioningProblem<inputAdapter_t> problem(&ia, &params);
191 #endif
192 
193  problem.solve();
194 }
195 
196 int main(int narg, char *arg[])
197 {
198  Tpetra::ScopeGuard tscope(&narg, &arg);
199  Teuchos::RCP<const Teuchos::Comm<int> > tcomm = Tpetra::getDefaultComm();
200 
201  int rank = tcomm->getRank();
202  int nParts = tcomm->getSize();
203  bool doRemap = false;
204  string filename = "USAir97";
205 
206  // Read run-time options.
207  Teuchos::CommandLineProcessor cmdp (false, false);
208  cmdp.setOption("file", &filename, "Name of the Matrix Market file to read");
209  cmdp.setOption("nparts", &nParts, "Number of parts.");
210  cmdp.setOption("remap", "no-remap", &doRemap, "Remap part numbers.");
211  cmdp.parse(narg, arg);
212 
213  meshCoordinatesTest(tcomm);
214 
215  testFromDataFile(tcomm, nParts, filename, doRemap);
216 
217  if (rank == 0)
218  serialTest(nParts, doRemap);
219 
220  if (rank == 0)
221  std::cout << "PASS" << std::endl;
222 }
void meshCoordinatesTest(const RCP< const Teuchos::Comm< int > > &comm)
Definition: rcbTest.cpp:163
int main(int narg, char *arg[])
double zscalar_t
static void getUIRandomData(unsigned int seed, zlno_t length, zscalar_t min, zscalar_t max, ArrayView< ArrayRCP< zscalar_t > > data)
Generate lists of random scalars.
Defines the PartitioningSolution class.
common code used by tests
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > myTypes_t
Definition: rcbTest.cpp:62
Defines the XpetraMultiVectorAdapter.
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
An adapter for Xpetra::MultiVector.
int zgno_t
PartitioningProblem sets up partitioning problems for the user.
void serialTest(int numParts, bool doRemap)
Definition: rcbTest.cpp:122
Tpetra::MultiVector< double, int, int > tMVector_t
RCP< tMVector_t > getUICoordinates()
#define nParts
Defines the PartitioningProblem class.
int testFromDataFile(RCP< const Teuchos::Comm< int > > &comm, int numParts, float imbalance, std::string fname, std::string pqParts, std::string pfname, int k, int migration_check_option, int migration_all_to_all_type, zscalar_t migration_imbalance_cut_off, int migration_processor_assignment_type, int migration_doMigration_type, bool test_boxes, bool rectilinear, int mj_premigration_option, int mj_premigration_coordinate_cutoff)
Defines the BasicVectorAdapter class.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
std::string testDataFilePath(".")