From 9a42b5484789ddfb00fb0d4421b7ab8881ef07ea Mon Sep 17 00:00:00 2001 From: Alexander Grayver Date: Mon, 5 Oct 2015 17:45:04 +0200 Subject: [PATCH] Make compute_intergrid_transfer_representation working for in parallel Update documentation --- include/deal.II/dofs/dof_tools.h | 8 +- source/dofs/dof_accessor_set.cc | 2 +- source/dofs/dof_tools_constraints.cc | 105 ++++++++++++++++++++------- 3 files changed, 84 insertions(+), 31 deletions(-) diff --git a/include/deal.II/dofs/dof_tools.h b/include/deal.II/dofs/dof_tools.h index 729dcb0420..61a0d54bca 100644 --- a/include/deal.II/dofs/dof_tools.h +++ b/include/deal.II/dofs/dof_tools.h @@ -740,8 +740,12 @@ namespace DoFTools * elements of the other vector components of the finite element fields on * the fine grid are not touched. * - * The output of this function is a compressed format that can be given to - * the @p reinit functions of the SparsityPattern and SparseMatrix classes. + * Triangulation of the fine grid can be distributed. When called in parallel, + * each process has to have a copy of the coarse grid. In this case, function + * returns transfer representation for a set of locally owned cells. + * + * The output of this function is a compressed format that can be used to + * construct corresponding sparse transfer matrix. */ template void diff --git a/source/dofs/dof_accessor_set.cc b/source/dofs/dof_accessor_set.cc index 6bdf85d14d..0260be236c 100644 --- a/source/dofs/dof_accessor_set.cc +++ b/source/dofs/dof_accessor_set.cc @@ -45,7 +45,7 @@ set_dof_values_by_interpolation (const Vector &local_values, OutputVector &values, const unsigned int fe_index) const { - if (!this->has_children()) + if (!this->has_children() && !this->is_artificial ()) { if ((dynamic_cast*> (this->dof_handler) diff --git a/source/dofs/dof_tools_constraints.cc b/source/dofs/dof_tools_constraints.cc index a07a13984d..81128949c8 100644 --- a/source/dofs/dof_tools_constraints.cc +++ b/source/dofs/dof_tools_constraints.cc @@ -34,6 +34,9 @@ #include #include +#ifdef DEAL_II_WITH_MPI +#include +#endif #include #include @@ -2500,7 +2503,11 @@ namespace DoFTools { unsigned int dofs_per_cell; std::vector parameter_dof_indices; +#ifdef DEAL_II_WITH_MPI + std::vector > global_parameter_representation; +#else std::vector > global_parameter_representation; +#endif }; } @@ -2568,7 +2575,6 @@ namespace DoFTools // vector to hold the representation of a single degree of freedom // on the coarse grid (for the selected fe) on the fine grid - const types::global_dof_index n_fine_dofs = weight_mapping.size(); copy_data.dofs_per_cell = coarse_fe.dofs_per_cell; copy_data.parameter_dof_indices.resize(copy_data.dofs_per_cell); @@ -2577,11 +2583,9 @@ namespace DoFTools // parameter grid cell cell->get_dof_indices (copy_data.parameter_dof_indices); - // reset the output array to a pristine state - copy_data.global_parameter_representation.clear (); - // loop over all dofs on this cell and check whether they are // interesting for us + unsigned int pos = 0; for (unsigned int local_dof=0; local_dof (n_fine_dofs)); + copy_data.global_parameter_representation[local_dof] = 0.; // distribute the representation of // @p{local_parameter_dof} on the parameter grid cell // @p{cell} to the global data space coarse_to_fine_grid_map[cell]-> set_dof_values_by_interpolation (parameter_dofs[local_parameter_dof], - copy_data.global_parameter_representation.back()); + copy_data.global_parameter_representation[pos]); + +#ifdef DEAL_II_WITH_MPI + copy_data.global_parameter_representation[pos].update_ghost_values (); +#endif + ++pos; } } @@ -2659,9 +2667,9 @@ namespace DoFTools weights[wi][wj] = copy_data.global_parameter_representation[pos](i); } } - else - Assert (copy_data.global_parameter_representation[pos](i) == 0, - ExcInternalError()); + //else + // Assert (copy_data.global_parameter_representation[pos](i) == 0, + // ExcInternalError()); ++pos; } @@ -2687,6 +2695,45 @@ namespace DoFTools Assembler::Scratch scratch; Assembler::CopyData copy_data; + unsigned int n_interesting_dofs = 0; + for (unsigned int local_dof=0; local_dof &tria = + dynamic_cast&> + (coarse_to_fine_grid_map.get_destination_grid().get_tria ()); + communicator = tria.get_communicator (); + } + catch (std::bad_cast &exp) + { + // Nothing bad happened: the user used serial Triangulation + } + + IndexSet locally_owned_dofs, locally_relevant_dofs; + DoFTools::extract_locally_owned_dofs + (coarse_to_fine_grid_map.get_destination_grid (), locally_owned_dofs); + DoFTools::extract_locally_relevant_dofs + (coarse_to_fine_grid_map.get_destination_grid (), locally_relevant_dofs); + + copy_data.global_parameter_representation[i].reinit + (locally_owned_dofs, locally_relevant_dofs, communicator); +#else + const types::global_dof_index n_fine_dofs = weight_mapping.size(); + copy_data.global_parameter_representation[i].resize (n_fine_dofs); +#endif + } + WorkStream::run(coarse_grid.begin_active(), coarse_grid.end(), std_cxx11::bind(&compute_intergrid_weights_3, @@ -2832,12 +2879,13 @@ namespace DoFTools for (typename dealii::DoFHandler::active_cell_iterator cell=fine_grid.begin_active(); cell!=fine_grid.end(); ++cell) - { - cell->get_dof_indices (local_dof_indices); - for (unsigned int i=0; iis_locally_owned ()) + { + cell->get_dof_indices (local_dof_indices); + for (unsigned int i=0; i::active_cell_iterator cell=fine_grid.begin_active(); cell != fine_grid.end(); ++cell) - { - cell->get_dof_indices (local_dof_indices); - for (unsigned int i=0; iis_locally_owned ()) + { + cell->get_dof_indices (local_dof_indices); + for (unsigned int i=0; i