From: kanschat Date: Mon, 18 Oct 2010 15:54:44 +0000 (+0000) Subject: no reason to be so picky about non-primitive elements X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=38c90dce56017cd841dd3b51f400d50633aab72c;p=dealii-svn.git no reason to be so picky about non-primitive elements git-svn-id: https://svn.dealii.org/trunk@22366 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 bedc72bb02..60fd39120d 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -1861,7 +1861,23 @@ class DoFTools /** * Make a constraint matrix for the * constraints that result from zero - * boundary values. This function is used + * boundary values. + * + * This function constrains all + * degrees of freedom on the + * boundary. Optionally, you can + * add a component mask, which + * restricts this functionality + * to a subset of an FESystem. + * + * For non-@ref GlossPrimitive "primitive" + * shape functions, any degree of freedom + * is affected that belongs to a + * shape function where at least + * one of its nonzero components + * is affected. + * + * This function is used * in step-36, for * example. * diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 574a8044f9..2a85824601 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -5770,27 +5770,7 @@ DoFTools::make_zero_boundary_constraints (const DH &dof, ++face_no) { const FiniteElement &fe = cell->get_fe(); - - // we can presently deal only with - // primitive elements for boundary - // values. make sure that all shape - // functions that are non-zero for - // the components we are interested - // in, are in fact primitive - for (unsigned int i=0; iget_fe().dofs_per_cell; ++i) - { - const std::vector &nonzero_component_array - = cell->get_fe().get_nonzero_components (i); - for (unsigned int c=0; cget_fe().is_primitive (i), - ExcMessage ("This function can only deal with requested boundary " - "values that correspond to primitive (scalar) base " - "elements")); - } - + typename DH::face_iterator face = cell->face(face_no); if (face->boundary_indicator () == 0) // face is of the right component @@ -5804,8 +5784,25 @@ DoFTools::make_zero_boundary_constraints (const DH &dof, // that match the component // signature. for (unsigned int i=0; i &nonzero_component_array + = cell->get_fe().get_nonzero_components (i); + bool nonzero = false; + for (unsigned int c=0; c