From: Martin Kronbichler Date: Fri, 10 Apr 2015 07:49:58 +0000 (+0200) Subject: Change arguments for SparsityTools::reorder_Cuthill_McKee and GridTools::get_face_con... X-Git-Tag: v8.3.0-rc1~299^2~5 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=2785bb952d5e78d70932c34f93b2d827419f72fa;p=dealii.git Change arguments for SparsityTools::reorder_Cuthill_McKee and GridTools::get_face_connectivity_of_cells to DynamicSparsityPattern. --- diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index e8429c9ea7..7483bb56f6 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -573,7 +573,17 @@ namespace GridTools template void get_face_connectivity_of_cells (const Triangulation &triangulation, - SparsityPattern &connectivity); + DynamicSparsityPattern &connectivity); + + /** + * As above, but filling a SparsityPattern object instead. + * + * @deprecated + */ + template + void + get_face_connectivity_of_cells (const Triangulation &triangulation, + SparsityPattern &connectivity) DEAL_II_DEPRECATED; /** * Use the METIS partitioner to generate a partitioning of the active cells diff --git a/include/deal.II/lac/sparsity_tools.h b/include/deal.II/lac/sparsity_tools.h index 09c05a4dde..1c5fb87248 100644 --- a/include/deal.II/lac/sparsity_tools.h +++ b/include/deal.II/lac/sparsity_tools.h @@ -20,6 +20,7 @@ #include #include #include +#include #include @@ -129,9 +130,19 @@ namespace SparsityTools * part of the algorithm that chooses starting indices. */ void - reorder_Cuthill_McKee (const SparsityPattern &sparsity, + reorder_Cuthill_McKee (const DynamicSparsityPattern &sparsity, + std::vector &new_indices, + const std::vector &starting_indices = std::vector()); + + /** + * As above, but taking a SparsityPattern object instead. + * + * @deprecated + */ + void + reorder_Cuthill_McKee (const SparsityPattern &sparsity, std::vector &new_indices, - const std::vector &starting_indices = std::vector()); + const std::vector &starting_indices = std::vector()) DEAL_II_DEPRECATED; #ifdef DEAL_II_WITH_MPI diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index f56143f2f2..7eb7533fe5 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -18,7 +18,7 @@ #include #include #include -#include +#include #include #include #include @@ -2489,7 +2489,7 @@ namespace parallel void Triangulation::setup_coarse_cell_to_p4est_tree_permutation () { - SparsityPattern cell_connectivity; + DynamicSparsityPattern cell_connectivity; GridTools::get_face_connectivity_of_cells (*this, cell_connectivity); coarse_cell_to_p4est_tree_permutation.resize (this->n_cells(0)); SparsityTools:: diff --git a/source/dofs/dof_renumbering.cc b/source/dofs/dof_renumbering.cc index f45484c5f1..e209c33cff 100644 --- a/source/dofs/dof_renumbering.cc +++ b/source/dofs/dof_renumbering.cc @@ -398,63 +398,48 @@ namespace DoFRenumbering constraints.close (); IndexSet locally_owned = dof_handler.locally_owned_dofs(); - SparsityPattern sparsity; // otherwise compute the Cuthill-McKee permutation - if (DH::dimension == 1) - { - sparsity.reinit (dof_handler.n_dofs(), - dof_handler.n_dofs(), - dof_handler.max_couplings_between_dofs()); - DoFTools::make_sparsity_pattern (dof_handler, sparsity, constraints); - sparsity.compress(); - } - else - { - DynamicSparsityPattern dsp (dof_handler.n_dofs(), - dof_handler.n_dofs(), - dof_handler.locally_owned_dofs()); - DoFTools::make_sparsity_pattern (dof_handler, dsp, constraints); - - // If the index set is not complete, need to get indices in local - // index space. - if (dof_handler.locally_owned_dofs().n_elements() != - dof_handler.locally_owned_dofs().size()) - { - // Create sparsity pattern from dsp by transferring its indices to - // processor-local index space and doing Cuthill-McKee there - std::vector row_lengths(locally_owned.n_elements()); - for (unsigned int i=0; i row_entries; - for (unsigned int i=0; i row_entries; + for (unsigned int i=0; i new_indices(sparsity.n_rows()); - SparsityTools::reorder_Cuthill_McKee (sparsity, new_indices, + std::vector new_indices(dsp.n_rows()); + SparsityTools::reorder_Cuthill_McKee (dsp, new_indices, starting_indices); if (reversed_numbering) diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index edd0ddfc9e..adf7da1a7d 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -1520,21 +1520,10 @@ next_cell: template void get_face_connectivity_of_cells (const Triangulation &triangulation, - SparsityPattern &cell_connectivity) + DynamicSparsityPattern &cell_connectivity) { - // as built in this function, we - // only consider face neighbors, - // which leads to a fixed number of - // entries per row (don't forget - // that each cell couples with - // itself, and that neighbors can - // be refined) cell_connectivity.reinit (triangulation.n_active_cells(), - triangulation.n_active_cells(), - GeometryInfo::faces_per_cell - * GeometryInfo::max_children_per_face - + - 1); + triangulation.n_active_cells()); // create a map pair -> SparsityPattern index // TODO: we are no longer using user_indices for this because we can get @@ -1548,22 +1537,13 @@ next_cell: cell != triangulation.end(); ++cell, ++index) indexmap[std::pair(cell->level(),cell->index())] = index; - // next loop over all cells and - // their neighbors to build the - // sparsity pattern. note that it's - // a bit hard to enter all the - // connections when a neighbor has - // children since we would need to - // find out which of its children - // is adjacent to the current - // cell. this problem can be - // omitted if we only do something - // if the neighbor has no children - // -- in that case it is either on - // the same or a coarser level than - // we are. in return, we have to - // add entries in both directions - // for both cells + // next loop over all cells and their neighbors to build the sparsity + // pattern. note that it's a bit hard to enter all the connections when a + // neighbor has children since we would need to find out which of its + // children is adjacent to the current cell. this problem can be omitted + // if we only do something if the neighbor has no children -- in that case + // it is either on the same or a coarser level than we are. in return, we + // have to add entries in both directions for both cells index = 0; for (typename dealii::internal::ActiveCellIterator >::type cell = triangulation.begin_active(); @@ -1581,12 +1561,22 @@ next_cell: cell_connectivity.add (other_index, index); } } + } + - // now compress the so-built connectivity pattern - cell_connectivity.compress (); + + template + void + get_face_connectivity_of_cells (const Triangulation &triangulation, + SparsityPattern &cell_connectivity) + { + DynamicSparsityPattern dsp; + get_face_connectivity_of_cells(triangulation, dsp); + cell_connectivity.copy_from(dsp); } + template void partition_triangulation (const unsigned int n_partitions, diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index 870714d4e3..63f0288e80 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -158,6 +158,11 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS Triangulation &, const bool); + template + void get_face_connectivity_of_cells + (const Triangulation &triangulation, + DynamicSparsityPattern &cell_connectivity); + template void get_face_connectivity_of_cells (const Triangulation &triangulation, diff --git a/source/lac/sparsity_tools.cc b/source/lac/sparsity_tools.cc index d28a31345e..1fba415d92 100644 --- a/source/lac/sparsity_tools.cc +++ b/source/lac/sparsity_tools.cc @@ -144,29 +144,20 @@ namespace SparsityTools * invalid_size_type indicates that a node has not been numbered yet), * pick a valid starting index among the as-yet unnumbered one. */ - SparsityPattern::size_type - find_unnumbered_starting_index (const SparsityPattern &sparsity, - const std::vector &new_indices) + DynamicSparsityPattern::size_type + find_unnumbered_starting_index (const DynamicSparsityPattern &sparsity, + const std::vector &new_indices) { { - SparsityPattern::size_type starting_point = numbers::invalid_size_type; - SparsityPattern::size_type min_coordination = sparsity.n_rows(); - for (SparsityPattern::size_type row=0; rowis_valid_entry() == false) - break; - // post-condition after loop: coordination, i.e. the number of - // entries in this row is now j-rowstart[row] - if (static_cast(j-sparsity.begin(row)) < - min_coordination) + if (sparsity.row_length(row) < min_coordination) { - min_coordination = j-sparsity.begin(row); + min_coordination = sparsity.row_length(row); starting_point = row; } } @@ -180,7 +171,7 @@ namespace SparsityTools // starting point, e.g. the first unnumbered one if (starting_point == numbers::invalid_size_type) { - for (SparsityPattern::size_type i=0; i &new_indices, - const std::vector &starting_indices) + reorder_Cuthill_McKee (const DynamicSparsityPattern &sparsity, + std::vector &new_indices, + const std::vector &starting_indices) { Assert (sparsity.n_rows() == sparsity.n_cols(), ExcDimensionMismatch (sparsity.n_rows(), sparsity.n_cols())); @@ -214,20 +206,20 @@ namespace SparsityTools // store the indices of the dofs renumbered in the last round. Default to // starting points - std::vector last_round_dofs (starting_indices); + std::vector last_round_dofs (starting_indices); // initialize the new_indices array with invalid values std::fill (new_indices.begin(), new_indices.end(), numbers::invalid_size_type); // delete disallowed elements - for (SparsityPattern::size_type i=0; i=sparsity.n_rows())) last_round_dofs[i] = numbers::invalid_size_type; std::remove_if (last_round_dofs.begin(), last_round_dofs.end(), - std::bind2nd(std::equal_to(), + std::bind2nd(std::equal_to(), numbers::invalid_size_type)); // now if no valid points remain: find dof with lowest coordination number @@ -237,81 +229,55 @@ namespace SparsityTools new_indices)); // store next free dof index - SparsityPattern::size_type next_free_number = 0; + DynamicSparsityPattern::size_type next_free_number = 0; // enumerate the first round dofs - for (SparsityPattern::size_type i=0; i!=last_round_dofs.size(); ++i) + for (DynamicSparsityPattern::size_type i=0; i!=last_round_dofs.size(); ++i) new_indices[last_round_dofs[i]] = next_free_number++; - // now do as many steps as needed to - // renumber all dofs + // now do as many steps as needed to renumber all dofs while (true) { - // store the indices of the dofs to be - // renumbered in the next round - std::vector next_round_dofs; - - // find all neighbors of the - // dofs numbered in the last - // round - for (SparsityPattern::size_type i=0; iis_valid_entry() == false) - break; - else - next_round_dofs.push_back (j->column()); + // store the indices of the dofs to be renumbered in the next round + std::vector next_round_dofs; + + // find all neighbors of the dofs numbered in the last round + for (DynamicSparsityPattern::size_type i=0; i::iterator end_sorted; + std::vector::iterator end_sorted; end_sorted = std::unique (next_round_dofs.begin(), next_round_dofs.end()); next_round_dofs.erase (end_sorted, next_round_dofs.end()); - // eliminate dofs which are - // already numbered + // eliminate dofs which are already numbered for (int s=next_round_dofs.size()-1; s>=0; --s) if (new_indices[next_round_dofs[s]] != numbers::invalid_size_type) next_round_dofs.erase (next_round_dofs.begin() + s); - // check whether there are - // any new dofs in the - // list. if there are none, - // then we have completely - // numbered the current - // component of the - // graph. check if there are - // as yet unnumbered - // components of the graph - // that we would then have to - // do next + // check whether there are any new dofs in the list. if there are + // none, then we have completely numbered the current component of the + // graph. check if there are as yet unnumbered components of the graph + // that we would then have to do next if (next_round_dofs.empty()) { if (std::find (new_indices.begin(), new_indices.end(), numbers::invalid_size_type) == new_indices.end()) - // no unnumbered - // indices, so we can - // leave now + // no unnumbered indices, so we can leave now break; - // otherwise find a valid - // starting point for the - // next component of the - // graph and continue - // with numbering that - // one. we only do so if - // no starting indices - // were provided by the - // user (see the - // documentation of this - // function) so produce - // an error if we got - // here and starting - // indices were given + // otherwise find a valid starting point for the next component of + // the graph and continue with numbering that one. we only do so + // if no starting indices were provided by the user (see the + // documentation of this function) so produce an error if we got + // here and starting indices were given Assert (starting_indices.empty(), ExcMessage ("The input graph appears to have more than one " "component, but as stated in the documentation " @@ -326,48 +292,36 @@ namespace SparsityTools - // store for each coordination - // number the dofs with these - // coordination number - std::multimap dofs_by_coordination; + // store for each coordination number the dofs with these coordination + // number + std::multimap dofs_by_coordination; - // find coordination number for - // each of these dofs - for (std::vector::iterator s=next_round_dofs.begin(); + // find coordination number for each of these dofs + for (std::vector::iterator s=next_round_dofs.begin(); s!=next_round_dofs.end(); ++s) { - SparsityPattern::size_type coordination = 0; - for (SparsityPattern::iterator j=sparsity.begin(*s); - jis_valid_entry() == false) - break; - else - ++coordination; - - // insert this dof at its - // coordination number - const std::pair new_entry (coordination, *s); + DynamicSparsityPattern::size_type coordination = 0; + for (DynamicSparsityPattern::row_iterator j=sparsity.row_begin(*s); + j new_entry (coordination, *s); dofs_by_coordination.insert (new_entry); } - // assign new DoF numbers to - // the elements of the present - // front: - std::multimap::iterator i; + // assign new DoF numbers to the elements of the present front: + std::multimap::iterator i; for (i = dofs_by_coordination.begin(); i!=dofs_by_coordination.end(); ++i) new_indices[i->second] = next_free_number++; - // after that: copy this round's - // dofs for the next round + // after that: copy this round's dofs for the next round last_round_dofs = next_round_dofs; } - // test for all indices - // numbered. this mostly tests - // whether the - // front-marching-algorithm (which - // Cuthill-McKee actually is) has - // reached all points. + // test for all indices numbered. this mostly tests whether the + // front-marching-algorithm (which Cuthill-McKee actually is) has reached + // all points. Assert ((std::find (new_indices.begin(), new_indices.end(), numbers::invalid_size_type) == new_indices.end()) @@ -376,6 +330,25 @@ namespace SparsityTools ExcInternalError()); } + + + void + reorder_Cuthill_McKee (const SparsityPattern &sparsity, + std::vector &new_indices, + const std::vector &starting_indices) + { + DynamicSparsityPattern dsp(sparsity.n_rows(), sparsity.n_cols()); + for (unsigned int row=0; rowis_valid_entry() ; ++it) + dsp.add(row, it->column()); + } + reorder_Cuthill_McKee(dsp, new_indices, starting_indices); + } + + + #ifdef DEAL_II_WITH_MPI template void distribute_sparsity_pattern(DSP_t &dsp, diff --git a/tests/bits/gerold_2.cc b/tests/bits/gerold_2.cc index b07676abc9..55cf7110cc 100644 --- a/tests/bits/gerold_2.cc +++ b/tests/bits/gerold_2.cc @@ -25,7 +25,7 @@ #include #include -#include +#include #include #include @@ -56,7 +56,7 @@ void LaplaceProblem::run () std::ifstream input_file(SOURCE_DIR "/gerold_1.inp"); grid_in.read_ucd(input_file); - SparsityPattern cell_connectivity; + DynamicSparsityPattern cell_connectivity; GridTools::get_face_connectivity_of_cells (triangulation, cell_connectivity); std::vector permutation(triangulation.n_active_cells()); diff --git a/tests/lac/sparsity_tools_01.cc b/tests/lac/sparsity_tools_01.cc index d99c32280b..a6f232ea19 100644 --- a/tests/lac/sparsity_tools_01.cc +++ b/tests/lac/sparsity_tools_01.cc @@ -21,7 +21,7 @@ #include "../tests.h" #include -#include +#include #include #include @@ -35,24 +35,20 @@ int main () deallog.depth_console(0); deallog.threshold_double(1.e-10); - SparsityPattern sp (4,4,4); + DynamicSparsityPattern dsp (4,4); for (unsigned int i=0; i<4; ++i) - sp.add (i,i); + dsp.add (i,i); - // create a graph with components - // 0,2 and 1,3 that are - // disconnected - sp.add (0,2); - sp.add (2,0); + // create a graph with components 0,2 and 1,3 that are disconnected + dsp.add (0,2); + dsp.add (2,0); - sp.add (1,3); - sp.add (3,1); - - sp.compress (); + dsp.add (1,3); + dsp.add (3,1); // now find permutation std::vector permutation(4); - SparsityTools::reorder_Cuthill_McKee (sp, permutation); + SparsityTools::reorder_Cuthill_McKee (dsp, permutation); for (unsigned int i=0; i