From 829516f61fed66fcf46780a5620ca6f1ba2f9b4f Mon Sep 17 00:00:00 2001 From: Nivesh Dommaraju Date: Sun, 26 Nov 2017 11:47:53 -0600 Subject: [PATCH] Add support for Zoltan graph partitioner Zoltan is used for graph partitioning. The following page for zoltan's guide has more information. http://www.cs.sandia.gov/zoltan/ug_html/ug_alg_phg.html Implementation: - 4 query functions, as need by zoltan, are added in sparsity_tools.cc - additional option has been added to "partition_triangulation" function to choose between "metis" or "zoltan" --- include/deal.II/grid/grid_tools.h | 40 ++-- include/deal.II/lac/sparsity_tools.h | 42 +++- source/base/mpi.cc | 10 + source/grid/grid_tools.cc | 15 +- source/grid/grid_tools.inst.in | 6 +- source/lac/sparsity_tools.cc | 331 +++++++++++++++++++++------ 6 files changed, 341 insertions(+), 103 deletions(-) diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index c0fb24880f..32eb210e2a 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -29,6 +29,7 @@ #include #include #include +#include DEAL_II_DISABLE_EXTRA_DIAGNOSTICS #include @@ -1439,25 +1440,30 @@ namespace GridTools DynamicSparsityPattern &connectivity); /** - * Use the METIS partitioner to generate a partitioning of the active cells - * making up the entire domain. After calling this function, the subdomain - * ids of all active cells will have values between zero and @p - * n_partitions-1. You can access the subdomain id of a cell by using + * Use graph partitioner to partition the active cells making up the entire domain. + * After calling this function, the subdomain ids of all active cells will have values + * between zero and @p n_partitions-1. You can access the subdomain id of a cell by using * cell-@>subdomain_id(). * - * This function will generate an error if METIS is not installed unless @p - * n_partitions is one. I.e., you can write a program so that it runs in the - * single-processor single-partition case without METIS installed, and only - * requires METIS when multiple partitions are required. + * Use the third argument to select between partitioning algorithms provided by METIS or ZOLTAN. + * METIS is the default partitioner. + * + * If deal.II was not installed with ZOLTAN or METIS, this function will generate an error + * when corresponding partition method is chosen, unless @p n_partitions is one. + * I.e., you can write a program so that it runs in the single-processor single-partition + * case without packages installed, and only requires them installed when + * multiple partitions are required. */ template void partition_triangulation (const unsigned int n_partitions, - Triangulation &triangulation); + Triangulation &triangulation, + const SparsityTools::Partitioner partitioner = SparsityTools::Partitioner::metis + ); /** * This function does the same as the previous one, i.e. it partitions a - * triangulation using METIS into a number of subdomains identified by the + * triangulation using a partitioning algorithm into a number of subdomains identified by the * cell-@>subdomain_id() flag. * * The difference to the previous function is the second argument, a @@ -1470,7 +1476,7 @@ namespace GridTools * triangulation. If the sparsity pattern contains an entry at position * $(i,j)$, then this means that cells $i$ and $j$ (in the order in which * they are traversed by active cell iterators) are to be considered - * connected; METIS will then try to partition the domain in such a way that + * connected; partitioning algorithm will then try to partition the domain in such a way that * (i) the subdomains are of roughly equal size, and (ii) a minimal number * of connections are broken. * @@ -1488,20 +1494,22 @@ namespace GridTools * default is to just consider face neighbors, not neighboring cells that * are connected by edges or vertices. While the latter couple when using * continuous finite elements, they are typically still closely connected in - * the neighborship graph, and METIS will not usually cut important - * connections in this case. However, if there are vertices in the mesh - * where many cells (many more than the common 4 or 6 in 2d and 3d, + * the neighborship graph, and partitioning algorithm + * will not usually cut important connections in this case. However, if there are + * vertices in the mesh where many cells (many more than the common 4 or 6 in 2d and 3d, * respectively) come together, then there will be a significant number of * cells that are connected across a vertex, but several degrees removed in * the connectivity graph built only using face neighbors. In a case like - * this, METIS may sometimes make bad decisions and you may want to build + * this, partitioning algorithm may sometimes make bad decisions and you may want to build * your own connectivity graph. */ template void partition_triangulation (const unsigned int n_partitions, const SparsityPattern &cell_connection_graph, - Triangulation &triangulation); + Triangulation &triangulation, + const SparsityTools::Partitioner partitioner = SparsityTools::Partitioner::metis + ); /** * Generates a partitioning of the active cells making up the entire domain diff --git a/include/deal.II/lac/sparsity_tools.h b/include/deal.II/lac/sparsity_tools.h index 2a108252d1..e76d6d4b8c 100644 --- a/include/deal.II/lac/sparsity_tools.h +++ b/include/deal.II/lac/sparsity_tools.h @@ -24,6 +24,7 @@ #include #include +#include #ifdef DEAL_II_WITH_MPI #include @@ -45,15 +46,33 @@ DEAL_II_NAMESPACE_OPEN */ namespace SparsityTools { + + /** + * Enumerator with options for partitioner + */ + enum class Partitioner + { + /** + * Use METIS partitioner. + */ + metis = 0, + /** + * Use ZOLTAN partitioner. + */ + zoltan + }; + + /** - * Use the METIS partitioner to generate a partitioning of the degrees of + * Use a partitioning algorithm to generate a partitioning of the degrees of * freedom represented by this sparsity pattern. In effect, we view this * sparsity pattern as a graph of connections between various degrees of * freedom, where each nonzero entry in the sparsity pattern corresponds to * an edge between two nodes in the connection graph. The goal is then to * decompose this graph into groups of nodes so that a minimal number of * edges are cut by the boundaries between node groups. This partitioning is - * done by METIS. Note that METIS can only partition symmetric sparsity + * done by METIS or ZOLTAN, depending upon which partitioner is chosen in the fourth argument. + * The default is METIS. Note that METIS and ZOLTAN can only partition symmetric sparsity * patterns, and that of course the sparsity pattern has to be square. We do * not check for symmetry of the sparsity pattern, since this is an * expensive operation, but rather leave this as the responsibility of @@ -63,10 +82,11 @@ namespace SparsityTools * zero and @p n_partitions-1 for each node (i.e. row or column of the * matrix). * - * This function will generate an error if METIS is not installed unless @p - * n_partitions is one. I.e., you can write a program so that it runs in the - * single-processor single-partition case without METIS installed, and only - * requires METIS when multiple partitions are required. + * If deal.II was not installed with packages ZOLTAN or METIS, this function will generate an error + * when corresponding partition method is chosen, unless @p n_partitions is one. + * I.e., you can write a program so that it runs in the single-processor single-partition + * case without the packages installed, and only requires them installed when + * multiple partitions are required. * * Note that the sparsity pattern itself is not changed by calling this * function. However, you will likely use the information generated by @@ -79,7 +99,9 @@ namespace SparsityTools */ void partition (const SparsityPattern &sparsity_pattern, const unsigned int n_partitions, - std::vector &partition_indices); + std::vector &partition_indices, + const Partitioner partitioner = Partitioner::metis + ); /** * For a given sparsity pattern, compute a re-enumeration of row/column @@ -242,6 +264,12 @@ namespace SparsityTools int, int, << "The array has size " << arg1 << " but should have size " << arg2); + /** + * Exception + */ + DeclExceptionMsg (ExcZOLTANNotInstalled, + "The function you called requires ZOLTAN, but you did not " + "configure deal.II with ZOLTAN or zoltan_cpp.h is not available."); } /** diff --git a/source/base/mpi.cc b/source/base/mpi.cc index f829c0d5fb..97ff462e2f 100644 --- a/source/base/mpi.cc +++ b/source/base/mpi.cc @@ -49,6 +49,10 @@ # include #endif +#ifdef DEAL_II_TRILINOS_WITH_ZOLTAN +# include +#endif + DEAL_II_NAMESPACE_OPEN @@ -351,6 +355,12 @@ namespace Utilities PetscPopSignalHandler(); #endif + //Initialize zoltan +#ifdef DEAL_II_TRILINOS_WITH_ZOLTAN + float version; + Zoltan_Initialize(argc, argv, &version); +#endif + #ifdef DEAL_II_WITH_P4EST //Initialize p4est and libsc components sc_init(MPI_COMM_WORLD, 0, 0, nullptr, SC_LP_SILENT); diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index cc6c40ff50..a7b9c47a0e 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -46,6 +46,7 @@ #include #include #include +#include #include #include @@ -2823,7 +2824,9 @@ next_cell: template void partition_triangulation (const unsigned int n_partitions, - Triangulation &triangulation) + Triangulation &triangulation, + const SparsityTools::Partitioner partitioner + ) { Assert ((dynamic_cast*> (&triangulation) @@ -2856,7 +2859,9 @@ next_cell: sp_cell_connectivity.copy_from(cell_connectivity); partition_triangulation (n_partitions, sp_cell_connectivity, - triangulation); + triangulation, + partitioner + ); } @@ -2865,7 +2870,9 @@ next_cell: void partition_triangulation (const unsigned int n_partitions, const SparsityPattern &cell_connection_graph, - Triangulation &triangulation) + Triangulation &triangulation, + const SparsityTools::Partitioner partitioner + ) { Assert ((dynamic_cast*> (&triangulation) @@ -2894,7 +2901,7 @@ next_cell: // of freedom (which is associated with a // cell) std::vector partition_indices (triangulation.n_active_cells()); - SparsityTools::partition (cell_connection_graph, n_partitions, partition_indices); + SparsityTools::partition (cell_connection_graph, n_partitions, partition_indices, partitioner); // finally loop over all cells and set the // subdomain ids diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index 7552215724..39cb48ebfe 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -267,12 +267,14 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS template void partition_triangulation (const unsigned int, - Triangulation &); + Triangulation &, + const SparsityTools::Partitioner); template void partition_triangulation (const unsigned int, const SparsityPattern &, - Triangulation &); + Triangulation &, + const SparsityTools::Partitioner); template void partition_triangulation_zorder (const unsigned int, diff --git a/source/lac/sparsity_tools.cc b/source/lac/sparsity_tools.cc index 1d8c3365f5..fd450bd659 100644 --- a/source/lac/sparsity_tools.cc +++ b/source/lac/sparsity_tools.cc @@ -22,6 +22,7 @@ #include #include #include +#include #ifdef DEAL_II_WITH_MPI #include @@ -39,15 +40,265 @@ extern "C" DEAL_II_ENABLE_EXTRA_DIAGNOSTICS #endif +#ifdef DEAL_II_TRILINOS_WITH_ZOLTAN +DEAL_II_DISABLE_EXTRA_DIAGNOSTICS +#include +DEAL_II_ENABLE_EXTRA_DIAGNOSTICS +#endif + +#include + DEAL_II_NAMESPACE_OPEN namespace SparsityTools { + namespace + { + void partition_metis (const SparsityPattern &sparsity_pattern, + const unsigned int n_partitions, + std::vector &partition_indices) + { + // Make sure that METIS is actually + // installed and detected +#ifndef DEAL_II_WITH_METIS + (void)sparsity_pattern; + AssertThrow (false, ExcMETISNotInstalled()); +#else + + // generate the data structures for + // METIS. Note that this is particularly + // simple, since METIS wants exactly our + // compressed row storage format. we only + // have to set up a few auxiliary arrays + idx_t + n = static_cast(sparsity_pattern.n_rows()), + ncon = 1, // number of balancing constraints (should be >0) + nparts = static_cast(n_partitions), // number of subdomains to create + dummy; // the numbers of edges cut by the + // resulting partition + + // We can not partition n items into more than n parts. METIS will + // generate non-sensical output (everything is owned by a single process) + // and complain with a message (but won't return an error code!): + // ***Cannot bisect a graph with 0 vertices! + // ***You are trying to partition a graph into too many parts! + nparts = std::min(n, nparts); + + // use default options for METIS + idx_t options[METIS_NOPTIONS]; + METIS_SetDefaultOptions (options); + + // one more nuisance: we have to copy our own data to arrays that store + // signed integers :-( + std::vector int_rowstart(1); + int_rowstart.reserve(sparsity_pattern.n_rows()+1); + std::vector int_colnums; + int_colnums.reserve(sparsity_pattern.n_nonzero_elements()); + for (SparsityPattern::size_type row=0; rowcolumn()); + int_rowstart.push_back(int_colnums.size()); + } + + std::vector int_partition_indices (sparsity_pattern.n_rows()); + + // Make use of METIS' error code. + int ierr; + + // Select which type of partitioning to create + + // Use recursive if the number of partitions is less than or equal to 8 + if (nparts <= 8) + ierr = METIS_PartGraphRecursive(&n, &ncon, int_rowstart.data(), int_colnums.data(), + nullptr, nullptr, nullptr, + &nparts,nullptr,nullptr,options, + &dummy,int_partition_indices.data()); + + // Otherwise use kway + else + ierr = METIS_PartGraphKway(&n, &ncon, int_rowstart.data(), int_colnums.data(), + nullptr, nullptr, nullptr, + &nparts,nullptr,nullptr,options, + &dummy,int_partition_indices.data()); + + // If metis returns normally, an error code METIS_OK=1 is returned from + // the above functions (see metish.h) + AssertThrow (ierr == 1, ExcMETISError (ierr)); + + // now copy back generated indices into the output array + std::copy (int_partition_indices.begin(), + int_partition_indices.end(), + partition_indices.begin()); +#endif + } + + +//Query functions unused if zoltan is not installed +#ifdef DEAL_II_TRILINOS_WITH_ZOLTAN + //Query functions for partition_zoltan + int get_number_of_objects(void *data, int *ierr) + { + SparsityPattern *graph = reinterpret_cast(data); + + *ierr = ZOLTAN_OK; + + return graph->n_rows(); + } + + + void get_object_list(void *data, int sizeGID, int sizeLID, + ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID, + int wgt_dim, float *obj_wgts, int *ierr) + { + SparsityPattern *graph = reinterpret_cast(data); + *ierr = ZOLTAN_OK; + + Assert( globalID != nullptr , ExcInternalError() ); + Assert( localID != nullptr , ExcInternalError() ); + + //set global degrees of freedom + auto n_dofs = graph->n_rows(); + + for (unsigned int i=0; i < n_dofs; i++) + { + globalID[i] = i; //global ids run from 0 to n_dofs-1 + localID[i] = i; //same as global ids + } + } + + + void get_num_edges_list(void *data, int sizeGID, int sizeLID, + int num_obj, + ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID, + int *numEdges, int *ierr) + { + SparsityPattern *graph = reinterpret_cast(data); + + *ierr = ZOLTAN_OK; + + Assert ( numEdges != nullptr , ExcInternalError() ); + + for (unsigned int i=0; irow_length(globalID[i])-1; + } + + + + void get_edge_list(void *data, int sizeGID, int sizeLID, + int num_obj, ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID, + int *num_edges, + ZOLTAN_ID_PTR nborGID, int *nborProc, + int wgt_dim, float *ewgts, int *ierr) + { + SparsityPattern *graph = reinterpret_cast(data); + *ierr = ZOLTAN_OK; + + ZOLTAN_ID_PTR nextNborGID = nborGID; + int *nextNborProc = nborProc; + + Assert( nextNborGID != nullptr , ExcInternalError() ); + Assert( nextNborProc != nullptr , ExcInternalError() ); + + //Loop through rows corresponding to indices in globalID implicitly + for (SparsityPattern::size_type i=0; i < static_cast(num_obj); ++i) + { + //Loop through each column to find neighbours + for ( SparsityPattern::iterator col=graph->begin(i); + col < graph->end(i); ++col ) + //Ignore diagonal entries. Not needed for partitioning. + if ( i != col->column() ) + { + *nextNborGID++ = col->column(); + *nextNborProc++ = 0; //All the vertices on processor 0 + } + } + } +#endif + + + void partition_zoltan (const SparsityPattern &sparsity_pattern, + const unsigned int n_partitions, + std::vector &partition_indices) + { + // Make sure that ZOLTAN is actually + // installed and detected +#ifndef DEAL_II_TRILINOS_WITH_ZOLTAN + (void)sparsity_pattern; + AssertThrow (false, ExcZOLTANNotInstalled()); +#else + + //MPI environment must have been initialized by this point. + std::unique_ptr zz = std_cxx14::make_unique(MPI_COMM_SELF); + + //General parameters + zz->Set_Param( "LB_METHOD", "GRAPH" ); //graph based partition method (LB-load balancing) + zz->Set_Param( "DEBUG_LEVEL", "0" ); //set level of debug info + zz->Set_Param( "NUM_LOCAL_PARTS", std::to_string(n_partitions) ); //set number of partitions + + //Need a non-const object equal to sparsity_pattern + SparsityPattern graph; + graph.copy_from(sparsity_pattern); + + //Set query functions + zz->Set_Num_Obj_Fn(get_number_of_objects, &graph); + zz->Set_Obj_List_Fn(get_object_list, &graph); + zz->Set_Num_Edges_Multi_Fn(get_num_edges_list, &graph); + zz->Set_Edge_List_Multi_Fn(get_edge_list, &graph); + + //Variables needed by partition function + int changes=0; + int num_gid_entries=1; + int num_lid_entries=1; + int num_import=0; + ZOLTAN_ID_PTR import_global_ids=nullptr; + ZOLTAN_ID_PTR import_local_ids=nullptr; + int *import_procs=nullptr; + int *import_to_part=nullptr; + int num_export=0; + ZOLTAN_ID_PTR export_global_ids=nullptr; + ZOLTAN_ID_PTR export_local_ids=nullptr; + int *export_procs=nullptr; + int *export_to_part=nullptr; + + //call partitioner + int rc = zz->LB_Partition (changes, + num_gid_entries, + num_lid_entries, + num_import, + import_global_ids, + import_local_ids, + import_procs, + import_to_part, + num_export, + export_global_ids, + export_local_ids, + export_procs, + export_to_part); + + //check for error code in partitioner + Assert( rc == ZOLTAN_OK , ExcInternalError() ); + + //By default, all indices belong to part 0. After zoltan partition + //some are migrated to different part ID, which is stored in export_to_part array. + std::fill(partition_indices.begin(), partition_indices.end(), 0); + + //copy from export_to_part to partition_indices, whose part_ids != 0. + for (int i=0; i < num_export; i++) + partition_indices[export_local_ids[i]] = export_to_part[i]; +#endif + } + } + + void partition (const SparsityPattern &sparsity_pattern, const unsigned int n_partitions, - std::vector &partition_indices) + std::vector &partition_indices, + const Partitioner partitioner + ) { Assert (sparsity_pattern.n_rows()==sparsity_pattern.n_cols(), ExcNotQuadratic()); @@ -66,80 +317,12 @@ namespace SparsityTools return; } - // Make sure that METIS is actually - // installed and detected -#ifndef DEAL_II_WITH_METIS - (void)sparsity_pattern; - AssertThrow (false, ExcMETISNotInstalled()); -#else - - // generate the data structures for - // METIS. Note that this is particularly - // simple, since METIS wants exactly our - // compressed row storage format. we only - // have to set up a few auxiliary arrays - idx_t - n = static_cast(sparsity_pattern.n_rows()), - ncon = 1, // number of balancing constraints (should be >0) - nparts = static_cast(n_partitions), // number of subdomains to create - dummy; // the numbers of edges cut by the - // resulting partition - - // We can not partition n items into more than n parts. METIS will - // generate non-sensical output (everything is owned by a single process) - // and complain with a message (but won't return an error code!): - // ***Cannot bisect a graph with 0 vertices! - // ***You are trying to partition a graph into too many parts! - nparts = std::min(n, nparts); - - // use default options for METIS - idx_t options[METIS_NOPTIONS]; - METIS_SetDefaultOptions (options); - - // one more nuisance: we have to copy our own data to arrays that store - // signed integers :-( - std::vector int_rowstart(1); - int_rowstart.reserve(sparsity_pattern.n_rows()+1); - std::vector int_colnums; - int_colnums.reserve(sparsity_pattern.n_nonzero_elements()); - for (SparsityPattern::size_type row=0; rowcolumn()); - int_rowstart.push_back(int_colnums.size()); - } - - std::vector int_partition_indices (sparsity_pattern.n_rows()); - - // Make use of METIS' error code. - int ierr; - - // Select which type of partitioning to create - - // Use recursive if the number of partitions is less than or equal to 8 - if (nparts <= 8) - ierr = METIS_PartGraphRecursive(&n, &ncon, int_rowstart.data(), int_colnums.data(), - nullptr, nullptr, nullptr, - &nparts,nullptr,nullptr,options, - &dummy,int_partition_indices.data()); - - // Otherwise use kway + if (partitioner == Partitioner::metis) + partition_metis(sparsity_pattern, n_partitions, partition_indices); + else if (partitioner == Partitioner::zoltan) + partition_zoltan(sparsity_pattern, n_partitions, partition_indices); else - ierr = METIS_PartGraphKway(&n, &ncon, int_rowstart.data(), int_colnums.data(), - nullptr, nullptr, nullptr, - &nparts,nullptr,nullptr,options, - &dummy,int_partition_indices.data()); - - // If metis returns normally, an error code METIS_OK=1 is returned from - // the above functions (see metish.h) - AssertThrow (ierr == 1, ExcMETISError (ierr)); - - // now copy back generated indices into the output array - std::copy (int_partition_indices.begin(), - int_partition_indices.end(), - partition_indices.begin()); -#endif + AssertThrow(false, ExcInternalError()); } -- 2.39.5