From 39e0e494e1915e016ba13f641b5764a0f24ec78b Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 26 Oct 2000 15:48:17 +0000 Subject: [PATCH] Speed up the compute_integrid_constraints function by some parallelization. git-svn-id: https://svn.dealii.org/trunk@3461 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_tools.h | 38 +++ deal.II/deal.II/source/dofs/dof_tools.cc | 328 ++++++++++++++--------- 2 files changed, 246 insertions(+), 120 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index 1786321975..5b227216b9 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -870,6 +870,44 @@ class DoFTools * Exception */ DeclException0 (ExcInvalidBoundaryIndicator); + + private: + /** + * This is a helper function that + * is used in the computation of + * integrid constraints. See the + * function for a thorough + * description of how it works. + */ + template + static void + compute_intergrid_weights (const DoFHandler &coarse_grid, + const unsigned int coarse_component, + const InterGridMap &coarse_to_fine_grid_map, + const vector > ¶meter_dofs, + const vector &weight_mapping, + FullMatrix &weights); + + /** + * This is a function that is + * called by the previous one and + * that operates on a range of + * cells only. It is used to + * split up the whole range of + * cells into chunks which are + * then worked on in parallel, if + * multithreading is available. + */ + template + static void + compute_intergrid_weights_1 (const DoFHandler &coarse_grid, + const unsigned int coarse_component, + const InterGridMap &coarse_to_fine_grid_map, + const vector > ¶meter_dofs, + const vector &weight_mapping, + FullMatrix &weights, + const typename DoFHandler::active_cell_iterator &begin, + const typename DoFHandler::active_cell_iterator &end); }; diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index bc5f965742..9c6aa0b71d 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -12,6 +12,8 @@ //---------------------------- dof_tools.cc --------------------------- +#include +#include #include #include #include @@ -1155,11 +1157,6 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa true); }; - // vector to hold the representation of - // a single degree of freedom on the - // coarse grid (for the selected fe) - // on the fine grid - Vector global_parameter_representation (n_fine_dofs); // store the weights with which a dof // on the parameter grid contributes @@ -1219,6 +1216,201 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa }; + // for each cell on the parameter grid: + // find out which degrees of freedom on the + // fine grid correspond in which way to + // the degrees of freedom on the parameter + // grid + // + // do this in a separate function + // to allow for multithreading + // there. see this function also if + // you want to read more + // information on the algorithm + // used. + compute_intergrid_weights (coarse_grid, coarse_component, + coarse_to_fine_grid_map, parameter_dofs, + weight_mapping, weights); + + + // ok, now we have all weights for each + // dof on the fine grid. if in debug + // mode lets see if everything went smooth, + // i.e. each dof has sum of weights one + // + // in other words this means that + // if the sum of all shape + // functions on the parameter grid + // is one (which is always the + // case), then the representation + // on the state grid should be as + // well (division of unity) + // + // if the parameter grid has more + // than one component, then the + // respective dofs of the other + // components have sum of weights + // zero, of course. we do not + // explicitely ask which component + // a dof belongs to, but this at + // least tests some errors +#ifdef DEBUG + for (unsigned int col=0; col1) && (sum==0)), ExcInternalError()); + }; +#endif + + + // now we know that the weights in + // each row constitute a + // constraint. enter this into the + // constraints object + // + // first task: for each parameter + // dof on the parameter grid, find + // a representant on the fine, + // global grid. this is possible + // since we use conforming finite + // element. we take this + // representant to be the first + // element in this row with weight + // identical to one. the + // representant will become an + // unconstrained degree of freedom, + // while all others will be + // constrained to this dof (and + // possibly others) + vector representants(weights.m(), -1); + for (unsigned int parameter_dof=0; parameter_dof(column)) + break; + Assert (global_dof < weight_mapping.size(), ExcInternalError()); + + // now enter the representants global + // index into our list + representants[parameter_dof] = global_dof; + }; + + + // note for people that want to + // optimize this function: the + // largest part of the computing + // time is spent in the following, + // rather innocent block of + // code. basically, it must be the + // ConstraintMatrix::add_entry call + // which takes the bulk of the + // time, but it is not known to the + // author how to make it faster... + vector > constraint_line; + for (unsigned int global_dof=0; global_dof(global_dof))) + // dof unconstrained + continue; + + // otherwise enter all constraints + constraints.add_line (global_dof); + + constraint_line.clear (); + for (unsigned int row=first_used_row; row +void +DoFTools::compute_intergrid_weights (const DoFHandler &coarse_grid, + const unsigned int coarse_component, + const InterGridMap &coarse_to_fine_grid_map, + const vector > ¶meter_dofs, + const vector &weight_mapping, + FullMatrix &weights) +{ + // simply distribute the range of + // cells to different threads + typedef typename DoFHandler::active_cell_iterator active_cell_iterator; + vector > + cell_intervals = Threads::split_range (coarse_grid.begin_active(), + coarse_grid.end(), + multithread_info.n_default_threads); + + Threads::ThreadManager thread_manager; + for (unsigned int i=0; i) + .collect_args (coarse_grid, coarse_component, + coarse_to_fine_grid_map, parameter_dofs, + weight_mapping, weights, + cell_intervals[i].first, + cell_intervals[i].second)); + + // wait for the threads to finish + thread_manager.wait (); +}; + + + +template +void +DoFTools::compute_intergrid_weights_1 (const DoFHandler &coarse_grid, + const unsigned int coarse_component, + const InterGridMap &coarse_to_fine_grid_map, + const vector > ¶meter_dofs, + const vector &weight_mapping, + FullMatrix &weights, + const typename DoFHandler::active_cell_iterator &begin, + const typename DoFHandler::active_cell_iterator &end) +{ + // aliases to the finite elements + // used by the dof handlers: + const FiniteElement &coarse_fe = coarse_grid.get_fe(); + // for each cell on the parameter grid: // find out which degrees of freedom on the // fine grid correspond in which way to @@ -1283,10 +1475,18 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa // both cells 1 and 2, but the // correct weight is nevertheless // only 1. + + // vector to hold the representation of + // a single degree of freedom on the + // coarse grid (for the selected fe) + // on the fine grid + const unsigned int n_fine_dofs = weight_mapping.size(); + Vector global_parameter_representation (n_fine_dofs); + typename DoFHandler::active_cell_iterator cell; vector parameter_dof_indices (coarse_fe.dofs_per_cell); - for (cell=coarse_grid.begin_active(); cell!=coarse_grid.end(); ++cell) + for (cell=begin; cell!=end; ++cell) { // get the global indices of the // parameter dofs on this parameter @@ -1358,118 +1558,6 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa ExcInternalError()); }; }; - - // ok, now we have all weights for each - // dof on the fine grid. if in debug - // mode lets see if everything went smooth, - // i.e. each dof has sum of weights one - // - // in other words this means that - // if the sum of all shape - // functions on the parameter grid - // is one (which is always the - // case), then the representation - // on the state grid should be as - // well (division of unity) - // - // if the parameter grid has more - // than one component, then the - // respective dofs of the other - // components have sum of weights - // zero, of course. we do not - // explicitely ask which component - // a dof belongs to, but this at - // least tests some errors -#ifdef DEBUG - for (unsigned int col=0; col1) && (sum==0)), ExcInternalError()); - }; -#endif - - - // now we know that the weights in - // each row constitute a - // constraint. enter this into the - // constraints object - // - // first task: for each parameter - // dof on the parameter grid, find - // a representant on the fine, - // global grid. this is possible - // since we use conforming finite - // element. we take this - // representant to be the first - // element in this row with weight - // identical to one. the - // representant will become an - // unconstrained degree of freedom, - // while all others will be - // constrained to this dof (and - // possibly others) - vector representants(weights.m(), -1); - for (unsigned int parameter_dof=0; parameter_dof(column)) - break; - Assert (global_dof < weight_mapping.size(), ExcInternalError()); - - // now enter the representants global - // index into our list - representants[parameter_dof] = global_dof; - }; - - for (unsigned int global_dof=0; global_dof(global_dof))) - // dof unconstrained - continue; - - // otherwise enter all constraints - constraints.add_line (global_dof); - for (unsigned int row=first_used_row; row &dof_handler, const Vector &cell_data, - Vector &dof_data) const; + Vector &dof_data); template void DoFTools::distribute_cell_to_dof_vector (const DoFHandler &dof_handler, const Vector &cell_data, - Vector &dof_data) const; + Vector &dof_data); template void DoFTools::extract_dofs(const DoFHandler& dof, -- 2.39.5