From f1ffe103f15b29ca3ec79dc6ae8d865b16be0f71 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Thu, 19 Sep 2013 20:39:37 +0000 Subject: [PATCH] Change graph_coloring to accept an user-defined function to get the conflict_indices. git-svn-id: https://svn.dealii.org/trunk@30846 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/base/graph_coloring.h | 121 ++++++++++-------- tests/base/graph_coloring_01.cc | 14 +- tests/base/graph_coloring_02.cc | 13 +- tests/base/graph_coloring_03.cc | 13 +- 4 files changed, 107 insertions(+), 54 deletions(-) diff --git a/deal.II/include/deal.II/base/graph_coloring.h b/deal.II/include/deal.II/base/graph_coloring.h index f21b407a90..3bf6c6abf5 100644 --- a/deal.II/include/deal.II/base/graph_coloring.h +++ b/deal.II/include/deal.II/base/graph_coloring.h @@ -28,7 +28,7 @@ DEAL_II_NAMESPACE_OPEN -// Currently graph_coloring only works only for continuous FE. +/// This namespace contains the functions necessary to color a graph. namespace graph_coloring { /** @@ -37,38 +37,38 @@ namespace graph_coloring { */ template std::vector > create_partitioning(Iterator const &begin, - typename identity::type const &end) + typename identity::type const &end, + std::vector (*get_conflict_indices) (Iterator &)) { std::vector > partitioning(1,std::vector (1,begin)); - // Create a map from dofs to cells - boost::unordered_map > dof_to_cell; + // Number of iterators. + unsigned int n_iterators(0); + // Create a map from conflict indices to iterators + boost::unordered_map > indices_to_iterators; for (Iterator it=begin; it!=end; ++it) { - // Need to check that it works with hp - const unsigned int dofs_per_cell(it->get_fe().dofs_per_cell); - std::vector local_dof_indices(dofs_per_cell); - it->get_dof_indices(local_dof_indices); - for (unsigned int i=0; i conflict_indices = (*get_conflict_indices)(it); + const unsigned int n_conflict_indices(conflict_indices.size()); + for (unsigned int i=0; i used_it; used_it.insert(begin); - while (!done) + while (used_it.size()!=n_iterators) { typename std::vector::iterator vector_it(partitioning.back().begin()); typename std::vector::iterator vector_end(partitioning.back().end()); std::vector new_zone; for (; vector_it!=vector_end; ++vector_it) { - const unsigned int dofs_per_cell((*vector_it)->get_fe().dofs_per_cell); - std::vector local_dof_indices(dofs_per_cell); - (*vector_it)->get_dof_indices(local_dof_indices); - for (unsigned int i=0; i conflict_indices = (*get_conflict_indices)(*vector_it); + const unsigned int n_conflict_indices(conflict_indices.size()); + for (unsigned int i=0; i iterator_vector(dof_to_cell[local_dof_indices[i]]); + std::vector iterator_vector(indices_to_iterators[conflict_indices[i]]); for (unsigned int j=0; j (1,it)); + break; + } } return partitioning; } + + /** - * This function uses DSATUR (Degree SATURation) to color the graph: + * This function uses DSATUR (Degree SATURation) to color one zone of the + * partition. DSATUR works as follows: * -# Arrange the vertices by decreasing order of degrees. * -# Color a vertex of maximal degree with color 1. * -# Choose a vertex with a maximal saturation degree. If there is equality, @@ -101,39 +110,40 @@ namespace graph_coloring { * -# If all the vertices are colored, stop. Otherwise, return to 3. */ template - std::vector > make_dsatur_coloring(std::vector const &partition) + std::vector > make_dsatur_coloring(std::vector &partition, + std::vector (*get_conflict_indices)(Iterator &)) { std::vector > partition_coloring; // Number of zones composing the partitioning. const unsigned int partition_size(partition.size()); std::vector sorted_vertices(partition_size); std::vector degrees(partition_size); - std::vector > local_dof_indices(partition_size); + std::vector > conflict_indices(partition_size); std::vector > graph(partition_size); - // Get the dofs associated to each cell. The dof_indices have to be sorted so + // Get the conflict indices associated to each iterator. The conflict_indices have to be sorted so // set_intersection can be used later. for (unsigned int i=0; iget_fe().dofs_per_cell); - partition[i]->get_dof_indices(local_dof_indices[i]); - std::sort(local_dof_indices[i].begin(),local_dof_indices[i].end()); + conflict_indices[i] = (*get_conflict_indices)(partition[i]); + std::sort(conflict_indices[i].begin(),conflict_indices[i].end()); } - // Compute the degree of each vertex of the graph (cell) using the - // intersection of dofs - std::vector dofs_intersection; + // Compute the degree of each vertex of the graph using the + // intersection of the conflict indices. + std::vector conflict_indices_intersection; std::vector::iterator intersection_it; for (unsigned int i=0; i::iterator it; it = std::max_element(colors_counter[i].begin(),colors_counter[i].end()); - unsigned int min_cells(-1); + unsigned int min_iterators(-1); unsigned int pos(0); // Find the color of coloring with the least number of colors among // the colors that have not been used yet. for (unsigned int k=0; k::iterator it; it = std::max_element(colors_counter[i].begin(),colors_counter[i].end()); - unsigned int min_cells(-1); + unsigned int min_iterators(-1); unsigned int pos(0); // Find the color of coloring with the least number of colors among // the colors that have not been used yet. for (unsigned int k=0; k std::vector > - make_graph_coloring(Iterator const &begin,typename identity::type const &end) + make_graph_coloring(Iterator const &begin,typename identity::type const &end, + std::vector (*get_conflict_indices)(Iterator &)) { // Create the partitioning. - std::vector > partitioning = create_partitioning(begin,end); + std::vector > partitioning = create_partitioning(begin,end, + get_conflict_indices); - // Color the cells within each partition. + // Color the iterators within each partition. const unsigned int partitioning_size(partitioning.size()); std::vector > > partition_coloring( partitioning_size); for (unsigned int i=0; i #include +template +std::vector get_conflict_indices_cfem( + typename DoFHandler::active_cell_iterator &it) +{ + std::vector local_dof_indices(it->get_fe().dofs_per_cell); + it->get_dof_indices(local_dof_indices); + + return local_dof_indices; +} + template void check() { @@ -45,8 +55,10 @@ void check() dof_handler.distribute_dofs(fe); // Create the coloring + typename DoFHandler::active_cell_iterator cell(dof_handler.begin_active()); std::vector::active_cell_iterator> > coloring( - graph_coloring::make_graph_coloring(dof_handler.begin_active(),dof_handler.end())); + graph_coloring::make_graph_coloring(cell,dof_handler.end(), + get_conflict_indices_cfem)); // Output the coloring for (unsigned int color=0; color #include +template +std::vector get_conflict_indices_cfem( + typename DoFHandler::active_cell_iterator &it) +{ + std::vector local_dof_indices(it->get_fe().dofs_per_cell); + it->get_dof_indices(local_dof_indices); + + return local_dof_indices; +} + template void check() { @@ -57,7 +67,8 @@ void check() // Create the coloring std::vector::active_cell_iterator> > coloring( - graph_coloring::make_graph_coloring(dof_handler.begin_active(),dof_handler.end())); + graph_coloring::make_graph_coloring(dof_handler.begin_active(),dof_handler.end(), + get_conflict_indices_cfem)); // Output the coloring for (unsigned int color=0; color #include +template +std::vector get_conflict_indices_cfem( + typename hp::DoFHandler::active_cell_iterator &it) +{ + std::vector local_dof_indices(it->get_fe().dofs_per_cell); + it->get_dof_indices(local_dof_indices); + + return local_dof_indices; +} + template void check() { @@ -59,7 +69,8 @@ void check() // Create the coloring std::vector::active_cell_iterator> > coloring( - graph_coloring::make_graph_coloring(dof_handler.begin_active(),dof_handler.end())); + graph_coloring::make_graph_coloring(dof_handler.begin_active(),dof_handler.end(), + get_conflict_indices_cfem)); // Output the coloring for (unsigned int color=0; color