From: Wolfgang Bangerth Date: Thu, 3 May 2001 15:21:19 +0000 (+0000) Subject: Grant access to internal parts of a function. X-Git-Tag: v8.0.0~19259 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f715f73e19631fe5a6dc7d705080ff4cbbc9d2c1;p=dealii.git Grant access to internal parts of a function. git-svn-id: https://svn.dealii.org/trunk@4531 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index f20b26ee45..8d5222cfea 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -807,6 +807,7 @@ class DoFTools * only need one object of this * type. */ +//TODO:[WB] document last argument template static void compute_intergrid_constraints (const DoFHandler &coarse_grid, @@ -814,7 +815,8 @@ class DoFTools const DoFHandler &fine_grid, const unsigned int fine_component, const InterGridMap &coarse_to_fine_grid_map, - ConstraintMatrix &constraints); + ConstraintMatrix &constraints, + std::vector > *transfer_representation = 0); /** * Create a mapping from degree diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index c6e3a3f77d..d41c0ea1a6 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -1170,7 +1170,8 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa const DoFHandler &fine_grid, const unsigned int fine_component, const InterGridMap &coarse_to_fine_grid_map, - ConstraintMatrix &constraints) + ConstraintMatrix &constraints, + std::vector > *transfer_representation) { // aliases to the finite elements // used by the dof handlers: @@ -1333,7 +1334,7 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa // the weights array was actually // of FullMatrix type. this // wasted huge amounts of memory, - // but was fast. nontheless, since + // but was fast. nonetheless, since // the memory consumption was // quadratic in the number of // degrees of freedom, this was not @@ -1439,7 +1440,49 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &coa }; #endif + // if the user wants to have a + // representation of the transfer + // matrix, the provide it + if (transfer_representation != 0) + { + transfer_representation->clear (); + transfer_representation->resize (weights.size()); + + const unsigned int n_global_parm_dofs + = std::count_if (weight_mapping.begin(), weight_mapping.end(), + std::bind2nd (std::not_equal_to (), -1)); + + // first construct the inverse + // mapping of weight_mapping + std::vector inverse_weight_mapping (n_global_parm_dofs, + DoFHandler::invalid_dof_index); + for (unsigned int i=0; i(-1)) + { + Assert (parameter_dof < n_global_parm_dofs, ExcInternalError()); + Assert (inverse_weight_mapping[parameter_dof] == DoFHandler::invalid_dof_index, + ExcInternalError()); + + inverse_weight_mapping[parameter_dof] = i; + }; + }; + + // next copy over weights array + // and replace respective + // numbers + for (unsigned int i=0; i::const_iterator j = weights[i].begin(); + for (; j!=weights[i].end(); ++j) + (*transfer_representation)[i][inverse_weight_mapping[j->first]] = j->second; + }; + }; + // now we know that the weights in // each row constitute a // constraint. enter this into the @@ -2076,7 +2119,8 @@ DoFTools::compute_intergrid_constraints (const DoFHandler &, const DoFHandler &, const unsigned int , const InterGridMap &, - ConstraintMatrix &); + ConstraintMatrix &, + std::vector > *);