Zoltan 2 Version 0.5
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_BasicCoordinateInput.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::BasicCoordinateInput<myTypes> inputAdapter_t;
00117 
00119   // Create input data.
00120 
00121   size_t localCount = 40;
00122   int dim = 3;
00123 
00124   scalar_t *coords = new scalar_t [dim * localCount];
00125 
00126   scalar_t *x = coords; 
00127   scalar_t *y = x + localCount; 
00128   scalar_t *z = y + localCount; 
00129 
00130   // Create coordinates that range from 0 to 10.0
00131 
00132   srand(rank);
00133   scalar_t scalingFactor = 10.0 / RAND_MAX;
00134 
00135   for (size_t i=0; i < localCount*dim; i++){
00136     coords[i] = scalar_t(rand()) * scalingFactor;
00137   }
00138 
00139   // Create global ids for the coordinates.
00140 
00141   globalId_t *globalIds = new globalId_t [localCount];
00142   globalId_t offset = rank * localCount;
00143 
00144   for (size_t i=0; i < localCount; i++)
00145     globalIds[i] = offset++;
00146    
00148   // Create parameters for an RCB problem
00149 
00150   double tolerance = 1.1;
00151 
00152   if (rank == 0)
00153     std::cout << "Imbalance tolerance is " << tolerance << std::endl;
00154 
00155   Teuchos::ParameterList params("test params");
00156   params.set("debug_level", "basic_status");
00157   params.set("debug_procs", "0");
00158   params.set("error_check_level", "debug_mode_assertions");
00159 
00160   params.set("compute_metrics", "true");
00161   params.set("algorithm", "rcb");
00162   params.set("imbalance_tolerance", tolerance );
00163   params.set("num_global_parts", nprocs);
00164 
00165   params.set("bisection_num_test_cuts", 1);
00166    
00169   // A simple problem with no weights.
00172 
00173   // Create a Zoltan2 input adapter for this geometry. TODO explain
00174 
00175   inputAdapter_t ia1(localCount, globalIds, x, y, z, 1, 1, 1);
00176 
00177   // Create a Zoltan2 partitioning problem
00178 
00179 #ifdef HAVE_ZOLTAN2_MPI                   
00180   Zoltan2::PartitioningProblem<inputAdapter_t> *problem1 = 
00181            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia1, &params, 
00182                                                             MPI_COMM_WORLD);
00183 #else
00184   Zoltan2::PartitioningProblem<inputAdapter_t> *problem1 =
00185            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia1, &params);
00186 #endif
00187    
00188   // Solve the problem
00189 
00190   problem1->solve();
00191    
00192   // Check the solution.
00193 
00194   if (rank == 0)
00195     problem1->printMetrics(cout);
00196 
00197   if (rank == 0){
00198     scalar_t imb = problem1->getImbalance();
00199     if (imb <= tolerance)
00200       std::cout << "pass: " << imb << std::endl;
00201     else
00202       std::cout << "fail: " << imb << std::endl;
00203     std::cout << std::endl;
00204   }
00205    
00208   // Try a problem with weights (1 dimension)
00211 
00212   scalar_t *weights = new scalar_t [localCount];
00213   for (size_t i=0; i < localCount; i++){
00214     weights[i] = 1.0 + scalar_t(rank) / scalar_t(nprocs);
00215   }
00216 
00217   // Create a Zoltan2 input adapter that includes weights.
00218 
00219   vector<const scalar_t *>coordVec(2);
00220   vector<int> coordStrides(2);
00221 
00222   coordVec[0] = x; coordStrides[0] = 1;
00223   coordVec[1] = y; coordStrides[1] = 1;
00224 
00225   vector<const scalar_t *>weightVec(1);
00226   vector<int> weightStrides(1);
00227 
00228   weightVec[0] = weights; weightStrides[0] = 1;
00229 
00230   inputAdapter_t ia2(
00231     localCount, globalIds,  
00232     coordVec, coordStrides, 
00233     weightVec, weightStrides);
00234 
00235   // Create a Zoltan2 partitioning problem
00236 
00237 #ifdef HAVE_ZOLTAN2_MPI                   
00238   Zoltan2::PartitioningProblem<inputAdapter_t> *problem2 = 
00239            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia2, &params,
00240                                                             MPI_COMM_WORLD);
00241 #else
00242   Zoltan2::PartitioningProblem<inputAdapter_t> *problem2 =
00243            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia2, &params);
00244 #endif
00245 
00246   // Solve the problem
00247 
00248   problem2->solve();
00249 
00250   // Check the solution.
00251 
00252   if (rank == 0)
00253     problem2->printMetrics(cout);
00254 
00255   if (rank == 0){
00256     scalar_t imb = problem2->getImbalance();
00257     if (imb <= tolerance)
00258       std::cout << "pass: " << imb << std::endl;
00259     else
00260       std::cout << "fail: " << imb << std::endl;
00261     std::cout << std::endl;
00262   }
00263 
00264   if (localCount > 0){
00265     delete [] weights;
00266     weights = NULL;
00267   }
00268 
00271   // Try a problem with multiple weights.
00274 
00275   // Add to the parameters the multicriteria objective.
00276 
00277   params.set("partitioning_objective", "multicriteria_minimize_total_weight");
00278 
00279   // Create the new weights.
00280 
00281   weights = new scalar_t [localCount*3];
00282   srand(rank);
00283 
00284   for (size_t i=0; i < localCount*3; i+=3){
00285     weights[i] = 1.0 + rank / nprocs;      // weight dimension 1
00286     weights[i+1] = rank<nprocs/2 ? 1 : 2;  // weight dimension 2
00287     weights[i+2] = rand()/RAND_MAX +.5;    // weight dimension 3
00288   }
00289 
00290   // Create a Zoltan2 input adapter with these weights.
00291 
00292   weightVec.resize(3);
00293   weightStrides.resize(3);
00294 
00295   weightVec[0] = weights;   weightStrides[0] = 3;
00296   weightVec[1] = weights+1; weightStrides[1] = 3;
00297   weightVec[2] = weights+2; weightStrides[2] = 3;
00298 
00299   inputAdapter_t ia3(
00300     localCount, globalIds,  
00301     coordVec, coordStrides, 
00302     weightVec, weightStrides);
00303 
00304   // Create a Zoltan2 partitioning problem.
00305 
00306 #ifdef HAVE_ZOLTAN2_MPI                   
00307   Zoltan2::PartitioningProblem<inputAdapter_t> *problem3 = 
00308            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia3, &params,
00309                                                             MPI_COMM_WORLD);
00310 #else
00311   Zoltan2::PartitioningProblem<inputAdapter_t> *problem3 =
00312            new Zoltan2::PartitioningProblem<inputAdapter_t>(&ia3, &params);
00313 #endif
00314 
00315   // Solve the problem
00316 
00317   problem3->solve();
00318 
00319   // Check the solution.
00320 
00321   if (rank == 0)
00322     problem3->printMetrics(cout);
00323 
00324   if (rank == 0){
00325     scalar_t imb = problem3->getImbalance();
00326     if (imb <= tolerance)
00327       std::cout << "pass: " << imb << std::endl;
00328     else
00329       std::cout << "fail: " << imb << std::endl;
00330     std::cout << std::endl;
00331   }
00332 
00334   // Try the other multicriteria objectives.
00335 
00336   bool dataHasChanged = false;    // default is true
00337 
00338   params.set("partitioning_objective", "multicriteria_minimize_maximum_weight");
00339   problem3->resetParameters(&params);
00340   problem3->solve(dataHasChanged);    
00341   if (rank == 0){
00342     problem3->printMetrics(cout);
00343     scalar_t imb = problem3->getImbalance();
00344     if (imb <= tolerance)
00345       std::cout << "pass: " << imb << std::endl;
00346     else
00347       std::cout << "fail: " << imb << std::endl;
00348     std::cout << std::endl;
00349   }
00350 
00351   params.set("partitioning_objective", "multicriteria_balance_total_maximum");
00352   problem3->resetParameters(&params);
00353   problem3->solve(dataHasChanged);    
00354   if (rank == 0){
00355     problem3->printMetrics(cout);
00356     scalar_t imb = problem3->getImbalance();
00357     if (imb <= tolerance)
00358       std::cout << "pass: " << imb << std::endl;
00359     else
00360       std::cout << "fail: " << imb << std::endl;
00361     std::cout << std::endl;
00362   }
00363 
00364   if (localCount > 0){
00365     delete [] weights;
00366     weights = NULL;
00367   }
00368 
00371   // Using part sizes, ask for some parts to be empty.
00374 
00375   // Change the number of parts to twice the number of processes to
00376   // ensure that we have more than one global part.
00377 
00378   params.set("num_global_parts", nprocs*2);
00379 
00380   // Using the initial problem that did not have any weights, reset
00381   // parameter list, and give it some part sizes.
00382 
00383   problem1->resetParameters(&params);
00384 
00385   zoltan2_partId_t *partIds = new zoltan2_partId_t [2];
00386   scalar_t *partSizes = new scalar_t [2];
00387 
00388   partIds[0] = rank*2;    partSizes[0] = 0;
00389   partIds[1] = rank*2+1;  partSizes[1] = 1;
00390 
00391   problem1->setPartSizes(2, partIds, partSizes);
00392 
00393   // Solve the problem.  The argument "dataHasChanged" indicates 
00394   // that we have not changed the input data, which allows the problem
00395   // so skip some work when re-solving. 
00396 
00397   dataHasChanged = false;
00398 
00399   problem1->solve(dataHasChanged);
00400 
00401   // Obtain the solution
00402 
00403   const Zoltan2::PartitioningSolution<inputAdapter_t> &solution4 =
00404     problem1->getSolution();
00405 
00406   // Check it.  Part sizes should all be odd.
00407 
00408   const zoltan2_partId_t *partAssignments = solution4.getPartList();
00409 
00410   int numInEmptyParts = 0;
00411   for (size_t i=0; i < localCount; i++){
00412     if (partAssignments[i] % 2 == 0)
00413       numInEmptyParts++;
00414   }
00415 
00416   if (rank == 0)
00417     std::cout << "Request that " << nprocs << " parts be empty." <<std::endl;
00418 
00419   // Check the solution.
00420 
00421   if (rank == 0)
00422     problem1->printMetrics(cout);
00423 
00424   if (rank == 0){
00425     scalar_t imb = problem1->getImbalance();
00426     if (imb <= tolerance)
00427       std::cout << "pass: " << imb << std::endl;
00428     else
00429       std::cout << "fail: " << imb << std::endl;
00430     std::cout << std::endl;
00431   }
00432 
00433   delete [] partIds;
00434   partIds = NULL;
00435   delete [] partSizes;
00436   partSizes = NULL;
00437 
00438   if (coords)
00439     delete [] coords;
00440 
00441   if (globalIds)
00442     delete [] globalIds;
00443 
00444   delete problem1;
00445   delete problem2;
00446   delete problem3;
00447 
00448 #ifdef HAVE_ZOLTAN2_MPI
00449   MPI_Finalize();
00450 #endif
00451 
00452   if (rank == 0)
00453     std::cout << "PASS" << std::endl;
00454 }
00455