From: wolf Date: Fri, 7 May 1999 12:31:53 +0000 (+0000) Subject: Make distribute_cell_to_dof_vector aware of finite elements made up of several compon... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cc78d90f73783a4122ba6d4cf97d752f3d6e0545;p=dealii-svn.git Make distribute_cell_to_dof_vector aware of finite elements made up of several components. git-svn-id: https://svn.dealii.org/trunk@1299 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_handler.h b/deal.II/deal.II/include/dofs/dof_handler.h index 6a1fd77c0d..f5c6a12d65 100644 --- a/deal.II/deal.II/include/dofs/dof_handler.h +++ b/deal.II/deal.II/include/dofs/dof_handler.h @@ -722,10 +722,27 @@ class DoFHandler * The output vector, being a data * vector on the grid, always consists * of elements of type #double#. + * + * In case the finite element used by + * this DoFHandler consists of more than + * one component, you should give which + * component in the output vector should + * be used to store the finite element + * field in; the default is zero (no other + * value is allowed if the finite element + * consists only of one component). All + * other components of the vector remain + * untouched, i.e. their contents are + * not changed. + * + * It is assumed that the output vector + * #dof_data# already has the right size, + * i.e. #n_dofs()# elements. */ template void distribute_cell_to_dof_vector (const Vector &cell_data, - Vector &dof_data) const; + Vector &dof_data, + const unsigned int component = 0) const; /** * Create a mapping from degree of freedom @@ -1425,7 +1442,16 @@ class DoFHandler DeclException2 (ExcInvalidMaskDimension, int, int, << "The dimension of the mask " << arg1 << " does not match" - << " the number of components in the finite element object."); + << " the number of components in the finite element object."); + /** + * Exception + */ + DeclException2 (ExcInvalidComponent, + int, int, + << "The component you gave (" << arg1 << ") " + << "is invalid with respect to the number " + << "of components in the finite element " + << "(" << arg2 << ")"); protected: diff --git a/deal.II/deal.II/source/dofs/dof_handler.cc b/deal.II/deal.II/source/dofs/dof_handler.cc index fdb33a5dff..295d96993f 100644 --- a/deal.II/deal.II/source/dofs/dof_handler.cc +++ b/deal.II/deal.II/source/dofs/dof_handler.cc @@ -2553,12 +2553,21 @@ unsigned int DoFHandler<3>::max_couplings_between_boundary_dofs () const { template template void DoFHandler::distribute_cell_to_dof_vector (const Vector &cell_data, - Vector &dof_data) const { + Vector &dof_data, + const unsigned int component) const { Assert (cell_data.size()==tria->n_active_cells(), ExcWrongSize (cell_data.size(), tria->n_active_cells())); - - // assign the right size to the output vector - dof_data.reinit (n_dofs()); + Assert (dof_data.size()==n_dofs(), ExcWrongSize (dof_data.size(), n_dofs())); + Assert (component < selected_fe->n_components, + ExcInvalidComponent(component, selected_fe->n_components)); + + // store a flag whether we should care + // about different components. this is + // just a simplification, we could ask + // for this at every single place + // equally well + const bool consider_components = (selected_fe->n_components != 1); + // count how often we have added a value // in the sum for each dof vector touch_count (n_dofs(), 0); @@ -2568,28 +2577,40 @@ void DoFHandler::distribute_cell_to_dof_vector (const Vector &cell_ unsigned int present_cell = 0; const unsigned int total_dofs = selected_fe->total_dofs; vector dof_indices (total_dofs); - + for (; cell!=endc; ++cell, ++present_cell) { cell->get_dof_indices (dof_indices); for (unsigned int i=0; isystem_to_component_index(i).first == component)) + { + // sum up contribution of the + // present_cell to this dof + dof_data(dof_indices[i]) += cell_data(present_cell); + // note that we added another + // summand + ++touch_count[dof_indices[i]]; + }; }; - + // compute the mean value on all the // dofs by dividing with the number // of summands. for (unsigned int i=0; i