geometric/example_rcb.cpp

This example shows how to use geometric partitioning methods with multivectors.

00001 //@HEADER
00002 // ************************************************************************
00003 //
00004 //               Isorropia: Partitioning and Load Balancing Package
00005 //                 Copyright (2006) 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 //
00025 // ************************************************************************
00026 //@HEADER
00027 
00028 //
00029 // Read in the file "simple.coords" and perform geometric partitioning
00030 // with Zoltan's RCB method.
00031 //
00032 // --f={filename} will read a different coordinate file
00033 // --v will print out the partitioning (small coordinate files only)
00034 // --rib will use Recursive Inertial Bisection instead of RCB
00035 //
00036 // The input file should be a
00037 // text file containing 1, 2 or 3-dimensional coordinates, one per line,
00038 // with white space separating coordinate values.
00039 //
00040 // There are a few other coordinate files in the test/geometric directory.
00041 
00042 #include <Isorropia_ConfigDefs.hpp>
00043 #include <Isorropia_Epetra.hpp>
00044 #include <Isorropia_EpetraPartitioner.hpp>
00045 #include <Isorropia_EpetraRedistributor.hpp>
00046 
00047 #define PRINTLIMIT 5000
00048 
00049 #ifdef HAVE_EPETRA
00050 #include <Epetra_Import.h>
00051 #ifdef HAVE_MPI
00052 #include <Epetra_MpiComm.h>
00053 #else
00054 #include <Epetra_SerialComm.h>
00055 #endif
00056 #include <Epetra_MultiVector.h>
00057 
00058 #include <Teuchos_CommandLineProcessor.hpp>
00059 #include <Teuchos_RCP.hpp>
00060 
00061 #include <string>
00062 
00063 #include "ispatest_lbeval_utils.hpp"
00064 #include "ispatest_epetra_utils.hpp"
00065 #ifdef _MSC_VER
00066 #include "winprocess.h"
00067 #endif
00068 
00069 int main(int argc, char** argv) {
00070 
00071   int fail = 0, dim=0;  
00072   int localProc = 0;
00073   int numProcs = 1;
00074 
00075 #ifdef HAVE_MPI
00076   MPI_Init(&argc, &argv);
00077   MPI_Comm_rank(MPI_COMM_WORLD, &localProc);
00078   MPI_Comm_size(MPI_COMM_WORLD, &numProcs);
00079   const Epetra_MpiComm Comm(MPI_COMM_WORLD);
00080 #else
00081   const Epetra_SerialComm Comm;
00082 #endif
00083 
00084   if (getenv("DEBUGME")){
00085     std::cerr << localProc << " gdb test_rcb.exe " << getpid() << std::endl;
00086     sleep(15);
00087   }
00088 
00089   // =============================================================
00090   // get command line options
00091   // =============================================================
00092 
00093   Teuchos::CommandLineProcessor clp(false,true);
00094 
00095   std::string *inputFile = new std::string("simple.coords");
00096   bool verbose = false;
00097   bool rib = false;
00098 
00099   clp.setOption( "f", inputFile, "Name of coordinate input file");
00100 
00101   clp.setOption( "v", "q", &verbose,
00102                 "Display coordinates and weights before and after partitioning.");
00103 
00104   clp.setOption( "rib", "rcb", &rib, "Run RIB instead of RCB");
00105 
00106   Teuchos::CommandLineProcessor::EParseCommandLineReturn parse_return =
00107     clp.parse(argc,argv);
00108 
00109   if( parse_return == Teuchos::CommandLineProcessor::PARSE_HELP_PRINTED){
00110 #ifdef HAVE_MPI
00111     MPI_Finalize();
00112 #endif
00113     return 0;
00114   }
00115   if( parse_return != Teuchos::CommandLineProcessor::PARSE_SUCCESSFUL ) {
00116 #ifdef HAVE_MPI
00117     MPI_Finalize();
00118 #endif
00119     return 1;
00120   }
00121 
00122   // =============================================================
00123   // Open file of coordinates and distribute them across processes
00124   // =============================================================
00125   
00126   Epetra_MultiVector *mv = ispatest::file2multivector(Comm, *inputFile);
00127 
00128   if (!mv || ((dim = mv->NumVectors()) < 1)){
00129     if (localProc == 0)
00130       std::cerr << "Invalid input file " << *inputFile << std::endl;
00131     exit(1);
00132   }
00133 
00134   int vsize = mv->GlobalLength();
00135 
00136   if (verbose){
00137     if (vsize < PRINTLIMIT){
00138       ispatest::printMultiVector(*mv, std::cout, "Coordinates",PRINTLIMIT);
00139     }
00140     else{
00141       if (localProc == 0){
00142         std::cerr << "--v requested, but input file is larger than " << PRINTLIMIT << " coordinates." << std::endl;
00143         std::cerr << "Partitioning will be performed but will not be displayed." << std::endl;
00144         verbose = false;
00145       }
00146     }
00147   }
00148 
00149   // =============================================================
00150   // Create weights for coordinates - there are three different
00151   // functions in src/utils/ispatest_epetra_utils.cpp that can
00152   // create weights for the coordinates.
00153   // =============================================================
00154 
00155   //Epetra_MultiVector *wgts = ispatest::makeWeights(mv->Map(), &ispatest::unitWeights);
00156   Epetra_MultiVector *wgts = ispatest::makeWeights(mv->Map(), &ispatest::alternateWeights);
00157   //Epetra_MultiVector *wgts = ispatest::makeWeights(mv->Map(), &ispatest::veeWeights);
00158 
00159   Epetra_Vector * &w1 = (*wgts)(0);
00160 
00161   if (!wgts || ((dim = wgts->NumVectors()) != 1)){
00162     if (localProc == 0)
00163       std::cout << "can't create weights" << std::endl;
00164     exit(1);
00165   }
00166 
00167   if (verbose){ 
00168     ispatest::printMultiVector(*wgts, std::cout, "Weights",PRINTLIMIT);
00169   }
00170 
00171   // =============================================================
00172   //  Create a parameter list for Zoltan
00173   // =============================================================
00174 
00175   Teuchos::ParameterList params;
00176 
00177   if (rib){
00178 
00179    // This parameter specifies that partitioning should be performed
00180    // by a simple linear partitioner in Isorropia, rather than calling
00181    // upon Zoltan.
00182 
00183    params.set("Partitioning Method", "RIB");
00184    // params.set("PARTITIONING_METHOD", "RIB"); // same as above
00185   }
00186   else{
00187 
00188    // params.set("Partitioning Method", "RCB"); // default 
00189 
00190     // To set low-level Zoltan parameters, create a sublist titled "ZOLTAN"
00191     // with Zoltan parameters in it.  See the Zoltan Users' Guide for
00192     // Zoltan parameters.  http://www.cs.sandia.gov/Zoltan
00193 
00194     //Teuchos::ParameterList &sublist = params.sublist("ZOLTAN");
00195     //sublist.set("DEBUG_LEVEL", "1"); // Zoltan will print out parameters
00196     //sublist.set("DEBUG_LEVEL", "5");   // proc 0 will trace Zoltan calls
00197     //sublist.set("DEBUG_MEMORY", "2");  // Zoltan will trace alloc & free
00198   }
00199   // =============================================================
00200   // Create a partitioner, by default this will perform the partitioning as well
00201   // =============================================================
00202 
00203   Teuchos::RCP<const Epetra_MultiVector> mv_rcp = Teuchos::rcp(mv);
00204   Teuchos::RCP<const Epetra_MultiVector> wgts_rcp = Teuchos::rcp(wgts);
00205 
00206   Teuchos::RCP<Isorropia::Epetra::Partitioner> part =
00207     Teuchos::rcp(new Isorropia::Epetra::Partitioner(mv_rcp, wgts_rcp, params));
00208 
00209   // Create a Redistributor based on the partitioning
00210 
00211   Isorropia::Epetra::Redistributor rd(part);
00212 
00213   // Redistribute the coordinates
00214 
00215   Teuchos::RCP<Epetra_MultiVector> new_mv = rd.redistribute(*mv);
00216 
00217   if (verbose){
00218     ispatest::printMultiVector(*new_mv, std::cout, "New Coordinates", PRINTLIMIT);
00219   }
00220 
00221   // Redistribute the weights
00222 
00223   Teuchos::RCP<Epetra_MultiVector> new_wgts = rd.redistribute(*wgts);
00224 
00225   Epetra_Vector * &new_w1 = (*new_wgts)(0);
00226 
00227   if (verbose){
00228     ispatest::printMultiVector(*new_wgts, std::cout, "New Weights", PRINTLIMIT);
00229   }
00230 
00231   // =============================================================
00232   // Compute balance both before and after partitioning
00233   // =============================================================
00234 
00235   double min1, min2, max1, max2, avg1, avg2;
00236   double goal = 1.0 / (double)numProcs;
00237 
00238   ispatest::compute_balance(*w1, goal, min1, max1, avg1);
00239 
00240   if (localProc == 0){
00241     std::cout << "Balance before partitioning: min " ;
00242     std::cout << min1 << " max " << max1 << " average " << avg1 << std::endl;
00243   }
00244 
00245   ispatest::compute_balance(*new_w1, goal, min2, max2, avg2);
00246 
00247   if (localProc == 0){
00248     std::cout << "Balance after partitioning:  min ";
00249     std::cout << min2 << " max " << max2 << " average " << avg2 << std::endl;
00250   }
00251 
00252 #ifdef HAVE_MPI
00253   MPI_Finalize();
00254 #endif
00255 
00256   return fail;
00257 }
00258 
00259 #endif   // HAVE_EPETRA