Zoltan2
rcb_C.cpp
Go to the documentation of this file.
00001 // @HEADER
00002 //
00003 // ***********************************************************************
00004 //
00005 //   Zoltan2: A package of combinatorial algorithms for scientific computing
00006 //                  Copyright 2012 Sandia Corporation
00007 //
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Karen Devine      (kddevin@sandia.gov)
00039 //                    Erik Boman        (egboman@sandia.gov)
00040 //                    Siva Rajamanickam (srajama@sandia.gov)
00041 //
00042 // ***********************************************************************
00043 //
00044 // @HEADER
00045 
00050 #include <Zoltan2_PartitioningSolution.hpp>
00051 #include <Zoltan2_PartitioningProblem.hpp>
00052 #include <Zoltan2_BasicVectorAdapter.hpp>
00053 #include <Zoltan2_InputTraits.hpp>
00054 #include <vector>
00055 #include <cstdlib>
00056 
00057 using namespace std;
00058 using std::vector;
00059 
00064 // Zoltan2 is templated.  What data types will we use for
00065 // scalars (coordinate values and weights), for local ids, and
00066 // for global ids?
00067 //
00068 // If Zoltan2 was compiled with explicit instantiation, we will
00069 // use the library's data types.  These macros are defined
00070 // in Zoltan2_config.h.
00071 
00072 #ifdef HAVE_ZOLTAN2_INST_FLOAT_INT_LONG
00073 typedef float scalar_t;
00074 typedef int localId_t;
00075 typedef long globalId_t;
00076 #else
00077   #ifdef HAVE_ZOLTAN2_INST_DOUBLE_INT_LONG
00078   typedef double scalar_t;
00079   typedef int localId_t;
00080   typedef long globalId_t;
00081   #else
00082     #ifdef HAVE_ZOLTAN2_INST_FLOAT_INT_INT
00083     typedef float scalar_t;
00084     typedef int localId_t;
00085     typedef int globalId_t;
00086     #else
00087       #ifdef HAVE_ZOLTAN2_INST_DOUBLE_INT_INT
00088       typedef double scalar_t;
00089       typedef int localId_t;
00090       typedef int globalId_t;
00091       #else
00092       typedef float scalar_t;
00093       typedef int localId_t;
00094       typedef int globalId_t;
00095       #endif
00096     #endif
00097   #endif
00098 #endif
00099 
00100 
00101 int main(int argc, char *argv[])
00102 {
00103 #ifdef HAVE_ZOLTAN2_MPI                   
00104   MPI_Init(&argc, &argv);
00105   int rank, nprocs;
00106   MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
00107   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
00108 #else
00109   int rank=0, nprocs=1;
00110 #endif
00111 
00112   // TODO explain
00113   typedef Zoltan2::BasicUserTypes<scalar_t, globalId_t, localId_t, globalId_t> myTypes;
00114 
00115   // TODO explain
00116   typedef Zoltan2::BasicVectorAdapter<myTypes> inputAdapter_t;
00117   typedef inputAdapter_t::part_t part_t;
00118 
00120   // Create input data.
00121 
00122   size_t localCount = 40;
00123   int dim = 3;
00124 
00125   scalar_t *coords = new scalar_t [dim * localCount];
00126 
00127   scalar_t *x = coords; 
00128   scalar_t *y = x + localCount; 
00129   scalar_t *z = y + localCount; 
00130 
00131   // Create coordinates that range from 0 to 10.0
00132 
00133   srand(rank);
00134   scalar_t scalingFactor = 10.0 / RAND_MAX;
00135 
00136   for (size_t i=0; i < localCount*dim; i++){
00137     coords[i] = scalar_t(rand()) * scalingFactor;
00138   }
00139 
00140   // Create global ids for the coordinates.
00141 
00142   globalId_t *globalIds = new globalId_t [localCount];
00143   globalId_t offset = rank * localCount;
00144 
00145   for (size_t i=0; i < localCount; i++)
00146     globalIds[i] = offset++;
00147    
00149   // Create parameters for an RCB problem
00150 
00151   double tolerance = 1.1;
00152 
00153   if (rank == 0)
00154     std::cout << "Imbalance tolerance is " << tolerance << std::endl;
00155 
00156   Teuchos::ParameterList params("test params");
00157   params.set("debug_level", "basic_status");
00158   params.set("debug_procs", "0");
00159   params.set("error_check_level", "debug_mode_assertions");
00160 
00161   params.set("compute_metrics", "true");
00162   params.set("algorithm", "rcb");
00163   params.set("imbalance_tolerance", tolerance );
00164   params.set("num_global_parts", nprocs);
00165 
00166   params.set("bisection_num_test_cuts", 1);
00167    
00170   // A simple problem with no weights.
00173 
00174   // Create a Zoltan2 input adapter for this geometry. TODO explain
00175 
00176   inputAdapter_t ia1(localCount, globalIds, x, y, z, 1, 1, 1);
00177 
00178   // Create a Zoltan2 partitioning problem
00179 
00180 #ifdef HAVE_ZOLTAN2_MPI                   
00181   Zoltan2::PartitioningProblem<inputAdapter_t> *problem1 = 
00182            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia1, &params, 
00183                                                             MPI_COMM_WORLD);
00184 #else
00185   Zoltan2::PartitioningProblem<inputAdapter_t> *problem1 =
00186            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia1, &params);
00187 #endif
00188    
00189   // Solve the problem
00190 
00191   problem1->solve();
00192    
00193   // Check the solution.
00194 
00195   if (rank == 0)
00196     problem1->printMetrics(cout);
00197 
00198   if (rank == 0){
00199     scalar_t imb = problem1->getWeightImbalance();
00200     if (imb <= tolerance)
00201       std::cout << "pass: " << imb << std::endl;
00202     else
00203       std::cout << "fail: " << imb << std::endl;
00204     std::cout << std::endl;
00205   }
00206    
00209   // Try a problem with weights 
00212 
00213   scalar_t *weights = new scalar_t [localCount];
00214   for (size_t i=0; i < localCount; i++){
00215     weights[i] = 1.0 + scalar_t(rank) / scalar_t(nprocs);
00216   }
00217 
00218   // Create a Zoltan2 input adapter that includes weights.
00219 
00220   vector<const scalar_t *>coordVec(2);
00221   vector<int> coordStrides(2);
00222 
00223   coordVec[0] = x; coordStrides[0] = 1;
00224   coordVec[1] = y; coordStrides[1] = 1;
00225 
00226   vector<const scalar_t *>weightVec(1);
00227   vector<int> weightStrides(1);
00228 
00229   weightVec[0] = weights; weightStrides[0] = 1;
00230 
00231   inputAdapter_t ia2(
00232     localCount, globalIds,  
00233     coordVec, coordStrides, 
00234     weightVec, weightStrides);
00235 
00236   // Create a Zoltan2 partitioning problem
00237 
00238 #ifdef HAVE_ZOLTAN2_MPI                   
00239   Zoltan2::PartitioningProblem<inputAdapter_t> *problem2 = 
00240            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia2, &params,
00241                                                             MPI_COMM_WORLD);
00242 #else
00243   Zoltan2::PartitioningProblem<inputAdapter_t> *problem2 =
00244            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia2, &params);
00245 #endif
00246 
00247   // Solve the problem
00248 
00249   problem2->solve();
00250 
00251   // Check the solution.
00252 
00253   if (rank == 0)
00254     problem2->printMetrics(cout);
00255 
00256   if (rank == 0){
00257     scalar_t imb = problem2->getWeightImbalance();
00258     if (imb <= tolerance)
00259       std::cout << "pass: " << imb << std::endl;
00260     else
00261       std::cout << "fail: " << imb << std::endl;
00262     std::cout << std::endl;
00263   }
00264 
00265   if (localCount > 0){
00266     delete [] weights;
00267     weights = NULL;
00268   }
00269 
00272   // Try a problem with multiple weights.
00275 
00276   // Add to the parameters the multicriteria objective.
00277 
00278   params.set("partitioning_objective", "multicriteria_minimize_total_weight");
00279 
00280   // Create the new weights.
00281 
00282   weights = new scalar_t [localCount*3];
00283   srand(rank);
00284 
00285   for (size_t i=0; i < localCount*3; i+=3){
00286     weights[i] = 1.0 + rank / nprocs;      // weight idx 1
00287     weights[i+1] = rank<nprocs/2 ? 1 : 2;  // weight idx 2
00288     weights[i+2] = rand()/RAND_MAX +.5;    // weight idx 3
00289   }
00290 
00291   // Create a Zoltan2 input adapter with these weights.
00292 
00293   weightVec.resize(3);
00294   weightStrides.resize(3);
00295 
00296   weightVec[0] = weights;   weightStrides[0] = 3;
00297   weightVec[1] = weights+1; weightStrides[1] = 3;
00298   weightVec[2] = weights+2; weightStrides[2] = 3;
00299 
00300   inputAdapter_t ia3(
00301     localCount, globalIds,  
00302     coordVec, coordStrides, 
00303     weightVec, weightStrides);
00304 
00305   // Create a Zoltan2 partitioning problem.
00306 
00307 #ifdef HAVE_ZOLTAN2_MPI                   
00308   Zoltan2::PartitioningProblem<inputAdapter_t> *problem3 = 
00309            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia3, &params,
00310                                                             MPI_COMM_WORLD);
00311 #else
00312   Zoltan2::PartitioningProblem<inputAdapter_t> *problem3 =
00313            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia3, &params);
00314 #endif
00315 
00316   // Solve the problem
00317 
00318   problem3->solve();
00319 
00320   // Check the solution.
00321 
00322   if (rank == 0)
00323     problem3->printMetrics(cout);
00324 
00325   if (rank == 0){
00326     scalar_t imb = problem3->getWeightImbalance();
00327     if (imb <= tolerance)
00328       std::cout << "pass: " << imb << std::endl;
00329     else
00330       std::cout << "fail: " << imb << std::endl;
00331     std::cout << std::endl;
00332   }
00333 
00335   // Try the other multicriteria objectives.
00336 
00337   bool dataHasChanged = false;    // default is true
00338 
00339   params.set("partitioning_objective", "multicriteria_minimize_maximum_weight");
00340   problem3->resetParameters(&params);
00341   problem3->solve(dataHasChanged);    
00342   if (rank == 0){
00343     problem3->printMetrics(cout);
00344     scalar_t imb = problem3->getWeightImbalance();
00345     if (imb <= tolerance)
00346       std::cout << "pass: " << imb << std::endl;
00347     else
00348       std::cout << "fail: " << imb << std::endl;
00349     std::cout << std::endl;
00350   }
00351 
00352   params.set("partitioning_objective", "multicriteria_balance_total_maximum");
00353   problem3->resetParameters(&params);
00354   problem3->solve(dataHasChanged);    
00355   if (rank == 0){
00356     problem3->printMetrics(cout);
00357     scalar_t imb = problem3->getWeightImbalance();
00358     if (imb <= tolerance)
00359       std::cout << "pass: " << imb << std::endl;
00360     else
00361       std::cout << "fail: " << imb << std::endl;
00362     std::cout << std::endl;
00363   }
00364 
00365   if (localCount > 0){
00366     delete [] weights;
00367     weights = NULL;
00368   }
00369 
00372   // Using part sizes, ask for some parts to be empty.
00375 
00376   // Change the number of parts to twice the number of processes to
00377   // ensure that we have more than one global part.
00378 
00379   params.set("num_global_parts", nprocs*2);
00380 
00381   // Using the initial problem that did not have any weights, reset
00382   // parameter list, and give it some part sizes.
00383 
00384   problem1->resetParameters(&params);
00385 
00386   part_t partIds[2];
00387   scalar_t partSizes[2];
00388 
00389   partIds[0] = rank*2;    partSizes[0] = 0;
00390   partIds[1] = rank*2+1;  partSizes[1] = 1;
00391 
00392   problem1->setPartSizes(2, partIds, partSizes);
00393 
00394   // Solve the problem.  The argument "dataHasChanged" indicates 
00395   // that we have not changed the input data, which allows the problem
00396   // so skip some work when re-solving. 
00397 
00398   dataHasChanged = false;
00399 
00400   problem1->solve(dataHasChanged);
00401 
00402   // Obtain the solution
00403 
00404   const Zoltan2::PartitioningSolution<inputAdapter_t> &solution4 =
00405     problem1->getSolution();
00406 
00407   // Check it.  Part sizes should all be odd.
00408 
00409   const part_t *partAssignments = solution4.getPartList();
00410 
00411   int numInEmptyParts = 0;
00412   for (size_t i=0; i < localCount; i++){
00413     if (partAssignments[i] % 2 == 0)
00414       numInEmptyParts++;
00415   }
00416 
00417   if (rank == 0)
00418     std::cout << "Request that " << nprocs << " parts be empty." <<std::endl;
00419 
00420   // Check the solution.
00421 
00422   if (rank == 0)
00423     problem1->printMetrics(cout);
00424 
00425   if (rank == 0){
00426     scalar_t imb = problem1->getWeightImbalance();
00427     if (imb <= tolerance)
00428       std::cout << "pass: " << imb << std::endl;
00429     else
00430       std::cout << "fail: " << imb << std::endl;
00431     std::cout << std::endl;
00432   }
00433 
00434   if (coords)
00435     delete [] coords;
00436 
00437   if (globalIds)
00438     delete [] globalIds;
00439 
00440   delete problem1;
00441   delete problem2;
00442   delete problem3;
00443 
00444 #ifdef HAVE_ZOLTAN2_MPI
00445   MPI_Finalize();
00446 #endif
00447 
00448   if (rank == 0)
00449     std::cout << "PASS" << std::endl;
00450 }
00451