matrix_1.cpp

This example demonstrates creating a balanced copy of Epetra_CrsGraph and Epetra_CrsMatrix objects using Isorropia::Epetra::createBalancedCopy functions.

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 //This file is a self-contained example of creating an Epetra_CrsGraph
00030 //and Epetra_CrsMatrix object, and using Isorropia's createBalancedCopy
00031 //function to rebalance them.
00032 //--------------------------------------------------------------------
00033 
00034 //Include Isorropia_Exception.hpp only because the helper functions at
00035 //the bottom of this file (which create the epetra objects) can
00036 //potentially throw exceptions.
00037 #include <Isorropia_Exception.hpp>
00038 #include <Isorropia_EpetraCostDescriber.hpp>
00039 
00040 //The Isorropia user-interface functions being demonstrated are declared
00041 //in Isorropia_Epetra.hpp.
00042 #include <Isorropia_Epetra.hpp>
00043 
00044 #include <Teuchos_ParameterList.hpp>
00045 
00046 #ifdef HAVE_MPI
00047 #include <mpi.h>
00048 #endif
00049 
00050 #ifdef HAVE_EPETRA
00051 #ifdef HAVE_MPI
00052 #include <Epetra_MpiComm.h>
00053 #else
00054 #include <Epetra_SerialComm.h>
00055 #endif
00056 #include <Epetra_Map.h>
00057 #include <Epetra_CrsMatrix.h>
00058 #endif
00059 
00060 #include "ispatest_lbeval_utils.hpp"
00061 
00062 //Declarations for helper-functions that create epetra objects. These
00063 //functions are implemented at the bottom of this file.
00064 #ifdef HAVE_EPETRA
00065 Teuchos::RCP<Epetra_CrsGraph>
00066   create_epetra_graph(int numProcs, int localProc);
00067 
00068 Teuchos::RCP<Epetra_CrsMatrix>
00069   create_epetra_matrix(int numProcs, int localProc);
00070 #endif
00071 
00074 int main(int argc, char** argv) {
00075 #if defined(HAVE_MPI) && defined(HAVE_EPETRA)
00076 
00077   int numProcs = 1;
00078   int localProc = 0;
00079 
00080   //first, set up our MPI environment...
00081   MPI_Init(&argc, &argv);
00082   MPI_Comm_rank(MPI_COMM_WORLD, &localProc);
00083   MPI_Comm_size(MPI_COMM_WORLD, &numProcs);
00084 
00085   //Create a Epetra_CrsGraph object. This graph will be the input to
00086   //the Isorropia rebalancing function...
00087 
00088   Teuchos::RCP<Epetra_CrsGraph> crsgraph;
00089   try {
00090     crsgraph = create_epetra_graph(numProcs, localProc);
00091   }
00092   catch(std::exception& exc) {
00093     std::cout << "matrix_1 example: create_epetra_graph threw exception '"
00094           << exc.what() << "' on proc " << localProc << std::endl;
00095     MPI_Finalize();
00096     return(-1);
00097   }
00098 
00099   //Now call Isorropia::createBalancedCopyto create a balanced
00100   //copy of crsgraph. 
00101 
00102   if (localProc == 0) {
00103     std::cout << "Hypergraph partitioning" << std::endl;
00104   }
00105 
00106   Teuchos::ParameterList paramlist;
00107   //No parameters. By default, Isorropia will use Zoltan hypergraph 
00108   //partitioning, treating the graph columns as hyperedges and the
00109   //graph rows as vertices.
00110 
00111   Epetra_CrsGraph *balanced_graph = NULL;
00112   try {
00113     balanced_graph =
00114       Isorropia::Epetra::createBalancedCopy(*crsgraph, paramlist);
00115 
00116   }
00117   catch(std::exception& exc) {
00118     std::cout << "matrix_1 example: Isorropia::createBalancedCopy threw "
00119          << "exception '" << exc.what() << "' on proc "
00120          << localProc << std::endl;
00121     MPI_Finalize();
00122     return(-1);
00123   }
00124 
00125   // Results
00126 
00127   Isorropia::Epetra::CostDescriber emptyCostObject;
00128   double goalWeight = 1.0 / (double)numProcs;
00129   double bal0, bal1, cutn0, cutn1, cutl0, cutl1;
00130 
00131   // Balance and cut quality before partitioning
00132 
00133   ispatest::compute_hypergraph_metrics(*crsgraph, emptyCostObject, goalWeight,
00134                      bal0, cutn0, cutl0);
00135 
00136   // Balance and cut quality after partitioning
00137 
00138   ispatest::compute_hypergraph_metrics(*balanced_graph, emptyCostObject, goalWeight,
00139                      bal1, cutn1, cutl1);
00140 
00141   if (localProc == 0){
00142     std::cout << "Before partitioning hypergraph: ";
00143     std::cout << "Balance " << bal0 << " cutN " << cutn0 << " cutL " << cutl0;
00144     std::cout << std::endl;
00145 
00146     std::cout << "After partitioning hypergraph:  ";
00147     std::cout << "Balance " << bal1 << " cutN " << cutn1 << " cutL " << cutl1;
00148     std::cout << std::endl;
00149     std::cout << std::endl;
00150   }
00151 
00152 
00153   //Next, do a similar exercise with a Epetra_CrsMatrix. Like the
00154   //Epetra_CrsGraph example above, we'll create a matrix to use as input,
00155   //and then have Isorropia::createBalancedCopy create a copy which is
00156   //balanced so that the number of nonzeros are equal on each processor.
00157 
00158   Teuchos::RCP<Epetra_CrsMatrix> crsmatrix;
00159   try {
00160     crsmatrix = create_epetra_matrix(numProcs, localProc);
00161   }
00162   catch(std::exception& exc) {
00163     std::cout << "matrix_1 example: create_epetra_matrix threw exception '"
00164           << exc.what() << "' on proc " << localProc << std::endl;
00165     MPI_Finalize();
00166     return(-1);
00167   }
00168 
00169 #ifdef HAVE_ISORROPIA_ZOLTAN
00170   //This time, we'll try graph partitioning..
00171   //Create a parameter sublist for Zoltan parameters.
00172   paramlist.set("PARTITIONING METHOD", "GRAPH");
00173 
00174   //We could also call ParMetis, if Zoltan was built with ParMetis.
00175   //Teuchos::ParameterList& sublist = paramlist.sublist("Zoltan");
00176   //sublist.set("GRAPH_PACKAGE", "PARMETIS");
00177   //sublist.set("PARMETIS_METHOD", "PARTKWAY");
00178 
00179 #else
00180   //This should never happen, since Zoltan is now required by Isorropia.
00181 #endif
00182 
00183 
00184   if (localProc == 0) {
00185     std::cout << "Specifying GRAPH partitioning." << std::endl;
00186   }
00187 
00188   Epetra_CrsMatrix *balanced_matrix;
00189   try {
00190     balanced_matrix =
00191       Isorropia::Epetra::createBalancedCopy(*crsmatrix, paramlist);
00192   }
00193   catch(std::exception& exc) {
00194     std::cout << "matrix_1 example: Isorropia::createBalancedCopy(matrix)"
00195         << " threw exception '" << exc.what() << "' on proc "
00196          << localProc << std::endl;
00197     MPI_Finalize();
00198     return(-1);
00199   }
00200   // Results
00201 
00202   double cutWgt0, cutWgt1;
00203   int numCuts0, numCuts1;
00204 
00205   // Balance and cut quality before partitioning
00206 
00207   ispatest::compute_graph_metrics(*crsmatrix, emptyCostObject, goalWeight,
00208                      bal0, numCuts0, cutWgt0, cutn0, cutl0);
00209 
00210   // Balance and cut quality after partitioning
00211 
00212   ispatest::compute_graph_metrics(*balanced_matrix, emptyCostObject, goalWeight,
00213                      bal1, numCuts1, cutWgt1, cutn1, cutl1);
00214 
00215   if (localProc == 0){
00216     std::cout << "Before partitioning graph: Number of cuts " << numCuts0 << " Cut weight " << cutWgt0 << std::endl;
00217     std::cout << "                     Balance " << bal0 << " cutN " << cutn0 << " cutL " << cutl0;
00218     std::cout << std::endl;
00219 
00220     std::cout << "After partitioning graph:  Number of cuts " << numCuts1 << " Cut weight " << cutWgt1 << std::endl;
00221     std::cout << "                     Balance " << bal1 << " cutN " << cutn1 << " cutL " << cutl1;
00222     std::cout << std::endl;
00223   }
00224 
00225   MPI_Finalize();
00226 
00227 #else
00228   std::cout << "matrix_1: must have both MPI and EPETRA. Make sure Trilinos "
00229     << "is configured with --enable-mpi and --enable-epetra." << std::endl;
00230 #endif
00231 
00232   return(0);
00233 }
00234 
00235 //Below are implementations of the helper-functions that create the
00236 //poorly-balanced epetra objects for use in the above example program.
00237 
00238 #if defined(HAVE_MPI) && defined(HAVE_EPETRA)
00239 
00240 Teuchos::RCP<Epetra_CrsMatrix>
00241   create_epetra_matrix(int numProcs, int localProc)
00242 {
00243   //create an Epetra_CrsMatrix with rows spread un-evenly over
00244   //processors.
00245   Epetra_MpiComm comm(MPI_COMM_WORLD);
00246   int local_num_rows = 200;
00247   int nnz_per_row = local_num_rows/4+1;
00248   int global_num_rows = numProcs*local_num_rows;
00249 
00250   int mid_proc = numProcs/2;
00251   bool num_procs_even = numProcs%2==0 ? true : false;
00252 
00253   int adjustment = local_num_rows/2;
00254 
00255   //adjust local_num_rows so that it's not equal on all procs.
00256   if (localProc < mid_proc) {
00257     local_num_rows -= adjustment;
00258   }
00259   else {
00260     local_num_rows += adjustment;
00261   }
00262 
00263   //if numProcs is not an even number, undo the local_num_rows adjustment
00264   //on one proc so that the total will still be correct.
00265   if (localProc == numProcs-1) {
00266     if (num_procs_even == false) {
00267       local_num_rows -= adjustment;
00268     }
00269   }
00270 
00271   //now we're ready to create a row-map.
00272   Epetra_Map rowmap(global_num_rows, local_num_rows, 0, comm);
00273 
00274   //create a matrix
00275   Teuchos::RCP<Epetra_CrsMatrix> matrix =
00276     Teuchos::rcp(new Epetra_CrsMatrix(Copy, rowmap, nnz_per_row));
00277 
00278   std::vector<int> indices(nnz_per_row);
00279   std::vector<double> coefs(nnz_per_row);
00280 
00281   int err = 0;
00282 
00283   for(int i=0; i<local_num_rows; ++i) {
00284     int global_row = rowmap.GID(i);
00285     int first_col = global_row - nnz_per_row/2;
00286 
00287     if (first_col < 0) {
00288       first_col = 0;
00289     }
00290     else if (first_col > (global_num_rows - nnz_per_row)) {
00291       first_col = global_num_rows - nnz_per_row;
00292     }
00293 
00294     for(int j=0; j<nnz_per_row; ++j) {
00295       indices[j] = first_col + j;
00296       coefs[j] = 1.0;
00297     }
00298 
00299     int err = matrix->InsertGlobalValues(global_row, nnz_per_row,
00300                                          &coefs[0], &indices[0]);
00301     if (err < 0) {
00302       err = matrix->ReplaceGlobalValues(global_row, nnz_per_row,
00303                                         &coefs[0], &indices[0]);
00304       if (err < 0) {
00305         throw Isorropia::Exception("create_epetra_matrix: error inserting matrix values.");
00306       }
00307     }
00308   }
00309 
00310   err = matrix->FillComplete();
00311   if (err != 0) {
00312     throw Isorropia::Exception("create_epetra_matrix: error in matrix.FillComplete()");
00313   }
00314 
00315   return(matrix);
00316 }
00317 
00318 Teuchos::RCP<Epetra_CrsGraph>
00319   create_epetra_graph(int numProcs, int localProc)
00320 {
00321   if (localProc == 0) {
00322     std::cout << " creating Epetra_CrsGraph with un-even distribution..."
00323             << std::endl;
00324   }
00325 
00326   //create an Epetra_CrsGraph with rows spread un-evenly over
00327   //processors.
00328   Epetra_MpiComm comm(MPI_COMM_WORLD);
00329   int local_num_rows = 800;
00330   int nnz_per_row = local_num_rows/4+1;
00331   int global_num_rows = numProcs*local_num_rows;
00332 
00333   int mid_proc = numProcs/2;
00334   bool num_procs_even = numProcs%2==0 ? true : false;
00335 
00336   int adjustment = local_num_rows/2;
00337 
00338   //adjust local_num_rows so that it's not equal on all procs.
00339   if (localProc < mid_proc) {
00340     local_num_rows -= adjustment;
00341   }
00342   else {
00343     local_num_rows += adjustment;
00344   }
00345 
00346   //if numProcs is not an even number, undo the local_num_rows adjustment
00347   //on one proc so that the total will still be correct.
00348   if (localProc == numProcs-1) {
00349     if (num_procs_even == false) {
00350       local_num_rows -= adjustment;
00351     }
00352   }
00353 
00354   //now we're ready to create a row-map.
00355   Epetra_Map rowmap(global_num_rows, local_num_rows, 0, comm);
00356 
00357   //create a graph
00358   Teuchos::RCP<Epetra_CrsGraph> graph =
00359     Teuchos::rcp(new Epetra_CrsGraph(Copy, rowmap, nnz_per_row));
00360 
00361   std::vector<int> indices(nnz_per_row);
00362   std::vector<double> coefs(nnz_per_row);
00363 
00364   int err = 0;
00365 
00366   for(int i=0; i<local_num_rows; ++i) {
00367     int global_row = rowmap.GID(i);
00368     int first_col = global_row - nnz_per_row/2;
00369 
00370     if (first_col < 0) {
00371       first_col = 0;
00372     }
00373     else if (first_col > (global_num_rows - nnz_per_row)) {
00374       first_col = global_num_rows - nnz_per_row;
00375     }
00376 
00377     for(int j=0; j<nnz_per_row; ++j) {
00378       indices[j] = first_col + j;
00379       coefs[j] = 1.0;
00380     }
00381 
00382     err = graph->InsertGlobalIndices(global_row, nnz_per_row,
00383                                          &indices[0]);
00384     if (err < 0) {
00385       throw Isorropia::Exception("create_epetra_graph: error inserting indices in graph");
00386     }
00387   }
00388 
00389   err = graph->FillComplete();
00390   if (err != 0) {
00391     throw Isorropia::Exception("create_epetra_graph: error in graph.FillComplete()");
00392   }
00393 
00394   return(graph);
00395 }
00396 
00397 #endif //HAVE_MPI && HAVE_EPETRA
00398