From 494435ca8a9ceb5ba30d206f0a6afaf2f38d1e19 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 4 May 2006 16:47:26 +0000 Subject: [PATCH] Implement make_hanging_node_constraints for hp in the case of all the same elements on all cells. git-svn-id: https://svn.dealii.org/trunk@13053 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_tools.h | 118 ++--- deal.II/deal.II/source/dofs/dof_tools.cc | 635 ++++++++++++++--------- 2 files changed, 451 insertions(+), 302 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index 1ed08d628f..c264a55a95 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -590,77 +590,63 @@ class DoFTools */ /** - * Make up the constraints which - * is result from the use of hanging - * nodes. The object into which these - * are inserted is later - * used to condensate the global - * system matrices and to prolong - * the solution vectors from the true + * Compute the constraints + * resulting from the presence of + * hanging nodes. The object into + * which these are inserted is + * later used to condense the + * global system matrix and right + * hand side, and to extend the + * solution vectors from the true * degrees of freedom also to the - * constraint nodes. + * constraint nodes. This + * function is explained in + * detail in the @ref step_6 + * "step-6" tutorial program and + * is used in almost all + * following programs as well. * - * Since this method does not make sense in - * one dimension, the function returns - * immediately. The object is not cleared - * before use, so you should make sure - * it containts only constraints you still - * want; otherwise call the @p clear - * function. + * This function does not clear + * the constraint matrix object + * before use, in order to allow + * adding constraints from + * different sources to the same + * object. You therefore need to + * make sure it contains only + * constraints you still want; + * otherwise call the + * ConstraintMatrix::clear() + * function. Likewise, this + * function does not close the + * object since you may want to + * enter other constraints later + * on yourself. * - * To condense a given sparsity pattern, - * use ConstraintMatrix@p ::condense. - * Before doing so, you need to close - * the constraint object, which must be - * done after all constraints are entered. - * This function does not close the object - * since you may want to enter other - * constraints later on yourself. - */ - static void - make_hanging_node_constraints (const DoFHandler<1> &dof_handler, - ConstraintMatrix &constraints); - - /** - * Declaration of same function - * for different space dimension. - */ - static void - make_hanging_node_constraints (const DoFHandler<2> &dof_handler, - ConstraintMatrix &constraints); - - /** - * Declaration of same function - * for different space dimension. - */ - static void - make_hanging_node_constraints (const DoFHandler<3> &dof_handler, - ConstraintMatrix &constraints); - - /** - * Declaration of same function - * for hp::DoFHandler - */ - static void - make_hanging_node_constraints (const hp::DoFHandler<1> &dof_handler, - ConstraintMatrix &constraints); - - /** - * Declaration of same function - * for different space dimension. - */ - static void - make_hanging_node_constraints (const hp::DoFHandler<2> &dof_handler, - ConstraintMatrix &constraints); - - /** - * Declaration of same function - * for different space dimension. + * In the hp-case, i.e. when the + * argument is of type + * hp::DoFHandler, we consider + * constraints due to different + * finite elements used on two + * sides of a face between cells + * as hanging nodes as well. In + * other words, for hp finite + * elements, this function + * computes all constraints due + * to differing mesh sizes (h) or + * polynomial degrees (p) between + * adjacent cells. + * + * The template argument (and by + * consequence the type of the + * first argument to this + * function) can be either a + * ::DoFHandler, hp::DoFHandler, + * or MGDoFHandler. */ + template static void - make_hanging_node_constraints (const hp::DoFHandler<3> &dof_handler, - ConstraintMatrix &constraints); - + make_hanging_node_constraints (const DH &dof_handler, + ConstraintMatrix &constraints); //@} /** diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 6e2683de6a..3f84e789e5 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -1648,275 +1648,424 @@ DoFTools::make_flux_sparsity_pattern ( -#if deal_II_dimension == 1 - -void DoFTools::make_hanging_node_constraints ( - const DoFHandler<1> &, - ConstraintMatrix &) +namespace internal { - // nothing to be done here -} + namespace DoFTools + { + template class int2type + {}; -void DoFTools::make_hanging_node_constraints ( - const hp::DoFHandler<1> &, - ConstraintMatrix &) -{ - // we may have to compute - // constraints for vertices - Assert (false, ExcNotImplemented()); -} +#if deal_II_dimension == 1 + void + do_make_hanging_node_constraints (const ::DoFHandler<1> &, + ConstraintMatrix &, + int2type<1>) + { + // nothing to do for regular + // dof handlers in 1d + } -#endif + void + do_make_hanging_node_constraints (const ::hp::DoFHandler<1> &/*dof_handler*/, + ConstraintMatrix &/*constraints*/, + int2type<1>) + { + // we may have to compute + // constraints for + // vertices. gotta think about + // that a bit more +//TODO[WB]: think about what to do here... + } +#endif #if deal_II_dimension == 2 - -void DoFTools::make_hanging_node_constraints ( - const DoFHandler<2> &dof_handler, - ConstraintMatrix &constraints) -{ - const unsigned int dim = 2; - - const FiniteElement &fe = dof_handler.get_fe(); + template + void + do_make_hanging_node_constraints (const DH &dof_handler, + ConstraintMatrix &constraints, + int2type<2>) + { + const unsigned int dim = 2; - // have space for the degrees of - // freedom on mother and child - // lines - const unsigned int n_dofs_on_mother = 2*fe.dofs_per_vertex + fe.dofs_per_line, - n_dofs_on_children = fe.dofs_per_vertex + 2*fe.dofs_per_line; - - std::vector dofs_on_mother(n_dofs_on_mother); - std::vector dofs_on_children(n_dofs_on_children); - - Assert(n_dofs_on_mother == fe.constraints().n(), - ExcDimensionMismatch(n_dofs_on_mother, - fe.constraints().n())); - Assert(n_dofs_on_children == fe.constraints().m(), - ExcDimensionMismatch(n_dofs_on_children, - fe.constraints().m())); - - // loop over all lines; only on - // lines there can be constraints. - // We do so by looping over all - // active cells and checking - // whether any of the faces are - // refined which can only be from - // the neighboring cell because - // this one is active. In that - // case, the face is subject to - // constraints - // - // note that even though we may - // visit a face twice if the - // neighboring cells are equally - // refined, we can only visit each - // face with hanging nodes once - DoFHandler::active_cell_iterator cell = dof_handler.begin_active(), + std::vector dofs_on_mother; + std::vector dofs_on_children; + + // loop over all lines; only on + // lines there can be constraints. + // We do so by looping over all + // active cells and checking + // whether any of the faces are + // refined which can only be from + // the neighboring cell because + // this one is active. In that + // case, the face is subject to + // constraints + // + // note that even though we may + // visit a face twice if the + // neighboring cells are equally + // refined, we can only visit each + // face with hanging nodes once + typename DH::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end(); - for (; cell!=endc; ++cell) - for (unsigned int face=0; face::faces_per_cell; ++face) - if (cell->face(face)->has_children()) - { - const DoFHandler::line_iterator line = cell->face(face); - - // fill the dofs indices. Use same - // enumeration scheme as in - // @p{FiniteElement::constraints()} - unsigned int next_index = 0; - for (unsigned int vertex=0; vertex<2; ++vertex) - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - dofs_on_mother[next_index++] = line->vertex_dof_index(vertex,dof); - for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) - dofs_on_mother[next_index++] = line->dof_index(dof); - Assert (next_index == dofs_on_mother.size(), - ExcInternalError()); - - next_index = 0; - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - dofs_on_children[next_index++] = line->child(0)->vertex_dof_index(1,dof); - for (unsigned int child=0; child<2; ++child) - for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) - dofs_on_children[next_index++] = line->child(child)->dof_index(dof); - Assert (next_index == dofs_on_children.size(), - ExcInternalError()); - - // for each row in the constraint - // matrix for this line: - for (unsigned int row=0; row!=dofs_on_children.size(); ++row) + for (; cell!=endc; ++cell) + for (unsigned int face=0; face::faces_per_cell; ++face) + if (cell->face(face)->has_children()) { - constraints.add_line (dofs_on_children[row]); - for (unsigned int i=0; i!=dofs_on_mother.size(); ++i) - constraints.add_entry (dofs_on_children[row], - dofs_on_mother[i], - fe.constraints()(row,i)); - } - } -} - + // so now we've found a + // face of an active + // cell that has + // children. that means + // that there are + // hanging nodes here. + + // in any case, faces + // can have at most two + // active fe indices, + // but here the face + // can have only one + // (namely the same as + // that from the cell + // we're sitting on), + // and each of the + // children can have + // only one as + // well. check this + Assert (cell->face(face)->n_active_fe_indices() == 1, + ExcInternalError()); + Assert (cell->face(face)->fe_index_is_active(cell->active_fe_index()) + == true, + ExcInternalError()); + for (unsigned int c=0; c::subfaces_per_face; ++c) + Assert (cell->face(face)->child(c)->n_active_fe_indices() == 1, + ExcInternalError()); -void DoFTools::make_hanging_node_constraints ( - const hp::DoFHandler<2> &/*dof_handler*/, - ConstraintMatrix &/*constraints*/) -{ -//TODO[?]: Implement (required for continuous elements) - Assert (false, ExcNotImplemented()); -} + // right now, all that + // is implemented is + // the case that both + // sides use the same + // fe + for (unsigned int c=0; c::subfaces_per_face; ++c) + Assert (cell->face(face)->child(c) + ->fe_index_is_active(cell->active_fe_index()) == true, + ExcNotImplemented()); + + // ok, start up the work + const FiniteElement &fe = cell->get_fe(); + const unsigned int fe_index = cell->active_fe_index(); + + const unsigned int + n_dofs_on_mother = 2*fe.dofs_per_vertex + fe.dofs_per_line, + n_dofs_on_children = fe.dofs_per_vertex + 2*fe.dofs_per_line; + dofs_on_mother.resize (n_dofs_on_mother); + dofs_on_children.resize (n_dofs_on_children); + + Assert(n_dofs_on_mother == fe.constraints().n(), + ExcDimensionMismatch(n_dofs_on_mother, + fe.constraints().n())); + Assert(n_dofs_on_children == fe.constraints().m(), + ExcDimensionMismatch(n_dofs_on_children, + fe.constraints().m())); + + const typename DH::line_iterator this_face = cell->face(face); + + // fill the dofs indices. Use same + // enumeration scheme as in + // @p{FiniteElement::constraints()} + unsigned int next_index = 0; + for (unsigned int vertex=0; vertex<2; ++vertex) + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) + dofs_on_mother[next_index++] = this_face->vertex_dof_index(vertex,dof, + fe_index); + for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) + dofs_on_mother[next_index++] = this_face->dof_index(dof, fe_index); + Assert (next_index == dofs_on_mother.size(), + ExcInternalError()); + + next_index = 0; + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) + dofs_on_children[next_index++] + = this_face->child(0)->vertex_dof_index(1,dof,fe_index); + for (unsigned int child=0; child<2; ++child) + for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) + dofs_on_children[next_index++] + = this_face->child(child)->dof_index(dof, fe_index); + Assert (next_index == dofs_on_children.size(), + ExcInternalError()); + + // for each row in the constraint + // matrix for this line: + for (unsigned int row=0; row!=dofs_on_children.size(); ++row) + { + constraints.add_line (dofs_on_children[row]); + for (unsigned int i=0; i!=dofs_on_mother.size(); ++i) + constraints.add_entry (dofs_on_children[row], + dofs_on_mother[i], + fe.constraints()(row,i)); + } + } + else + { + // this face has no + // children, but it + // could still be that + // it is shared by two + // cells that use a + // different fe index + Assert (cell->face(face)->n_active_fe_indices() == 1, + ExcNotImplemented()); + Assert (cell->face(face) + ->fe_index_is_active(cell->active_fe_index()) == true, + ExcInternalError()); + } + } #endif #if deal_II_dimension == 3 - -void DoFTools::make_hanging_node_constraints ( - const DoFHandler<3> &dof_handler, - ConstraintMatrix &constraints) -{ - const unsigned int dim = 3; - - const FiniteElement &fe = dof_handler.get_fe(); + template + void + do_make_hanging_node_constraints (const DH &dof_handler, + ConstraintMatrix &constraints, + int2type<3>) + { + const unsigned int dim = 3; - // have space for the degrees of - // freedom on mother and child - // lines - const unsigned int - n_dofs_on_mother = (4*fe.dofs_per_vertex+ - 4*fe.dofs_per_line+ - fe.dofs_per_quad), - n_dofs_on_children = (5*fe.dofs_per_vertex+ - 12*fe.dofs_per_line+ - 4*fe.dofs_per_quad); - - std::vector dofs_on_mother(n_dofs_on_mother); - std::vector dofs_on_children(n_dofs_on_children); - - Assert(n_dofs_on_mother == fe.constraints().n(), - ExcDimensionMismatch(n_dofs_on_mother, - fe.constraints().n())); - Assert(n_dofs_on_children == fe.constraints().m(), - ExcDimensionMismatch(n_dofs_on_children, - fe.constraints().m())); - - // loop over all lines; only on - // lines there can be constraints. - // We do so by looping over all - // active cells and checking - // whether any of the faces are - // refined which can only be from - // the neighboring cell because - // this one is active. In that - // case, the face is subject to - // constraints - // - // note that even though we may - // visit a face twice if the - // neighboring cells are equally - // refined, we can only visit each - // face with hanging nodes once - DoFHandler::active_cell_iterator cell = dof_handler.begin_active(), - endc = dof_handler.end(); - for (; cell!=endc; ++cell) - for (unsigned int f=0; f::faces_per_cell; ++f) - if (cell->face(f)->has_children()) - { - const DoFHandler::face_iterator face = cell->face(f); - - // fill the dofs indices. Use same - // enumeration scheme as in - // @p{FiniteElement::constraints()} - unsigned int next_index = 0; - for (unsigned int vertex=0; vertex<4; ++vertex) - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - dofs_on_mother[next_index++] = face->vertex_dof_index(vertex,dof); - for (unsigned int line=0; line<4; ++line) - for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) - dofs_on_mother[next_index++] = face->line(line)->dof_index(dof); - for (unsigned int dof=0; dof!=fe.dofs_per_quad; ++dof) - dofs_on_mother[next_index++] = face->dof_index(dof); - Assert (next_index == dofs_on_mother.size(), - ExcInternalError()); + std::vector dofs_on_mother; + std::vector dofs_on_children; + + // loop over all quads; only on + // quads there can be constraints. + // We do so by looping over all + // active cells and checking + // whether any of the faces are + // refined which can only be from + // the neighboring cell because + // this one is active. In that + // case, the face is subject to + // constraints + // + // note that even though we may + // visit a face twice if the + // neighboring cells are equally + // refined, we can only visit each + // face with hanging nodes once + typename DH::active_cell_iterator cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell!=endc; ++cell) + for (unsigned int face=0; face::faces_per_cell; ++face) + if (cell->face(face)->has_children()) + { + // so now we've found a + // face of an active + // cell that has + // children. that means + // that there are + // hanging nodes here. + + // in any case, faces + // can have at most two + // active fe indices, + // but here the face + // can have only one + // (namely the same as + // that from the cell + // we're sitting on), + // and each of the + // children can have + // only one as + // well. check this + Assert (cell->face(face)->n_active_fe_indices() == 1, + ExcInternalError()); + Assert (cell->face(face)->fe_index_is_active(cell->active_fe_index()) + == true, + ExcInternalError()); + for (unsigned int c=0; c::subfaces_per_face; ++c) + Assert (cell->face(face)->child(c)->n_active_fe_indices() == 1, + ExcInternalError()); + + // right now, all that + // is implemented is + // the case that both + // sides use the same + // fe, and not only + // that but also that + // all lines bounding + // this face and the + // children have the + // same fe + for (unsigned int c=0; c::subfaces_per_face; ++c) + { + Assert (cell->face(face)->child(c) + ->fe_index_is_active(cell->active_fe_index()) == true, + ExcNotImplemented()); + for (unsigned int e=0; e<4; ++e) + { + Assert (cell->face(face)->child(c)->line(e) + ->n_active_fe_indices() == 1, + ExcNotImplemented()); + Assert (cell->face(face)->child(c)->line(e) + ->fe_index_is_active(cell->active_fe_index()) == true, + ExcNotImplemented()); + } + } + for (unsigned int e=0; e<4; ++e) + { + Assert (cell->face(face)->line(e) + ->n_active_fe_indices() == 1, + ExcNotImplemented()); + Assert (cell->face(face)->line(e) + ->fe_index_is_active(cell->active_fe_index()) == true, + ExcNotImplemented()); + } + + // ok, start up the work + const FiniteElement &fe = cell->get_fe(); + const unsigned int fe_index = cell->active_fe_index(); + + const unsigned int + n_dofs_on_mother = (4*fe.dofs_per_vertex+ + 4*fe.dofs_per_line+ + fe.dofs_per_quad), + n_dofs_on_children = (5*fe.dofs_per_vertex+ + 12*fe.dofs_per_line+ + 4*fe.dofs_per_quad); + + dofs_on_mother.resize (n_dofs_on_mother); + dofs_on_children.resize (n_dofs_on_children); + + Assert(n_dofs_on_mother == fe.constraints().n(), + ExcDimensionMismatch(n_dofs_on_mother, + fe.constraints().n())); + Assert(n_dofs_on_children == fe.constraints().m(), + ExcDimensionMismatch(n_dofs_on_children, + fe.constraints().m())); + + const typename DH::face_iterator this_face = cell->face(face); - next_index = 0; - - // assert some consistency - // assumptions - Assert ((face->child(0)->vertex_index(3) == - face->child(1)->vertex_index(2)) && - (face->child(0)->vertex_index(3) == - face->child(2)->vertex_index(1)) && - (face->child(0)->vertex_index(3) == - face->child(3)->vertex_index(0)), - ExcInternalError()); - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - dofs_on_children[next_index++] - = face->child(0)->vertex_dof_index(3,dof); + // fill the dofs indices. Use same + // enumeration scheme as in + // @p{FiniteElement::constraints()} + unsigned int next_index = 0; + for (unsigned int vertex=0; vertex<4; ++vertex) + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) + dofs_on_mother[next_index++] = this_face->vertex_dof_index(vertex,dof, + fe_index); + for (unsigned int line=0; line<4; ++line) + for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) + dofs_on_mother[next_index++] + = this_face->line(line)->dof_index(dof, fe_index); + for (unsigned int dof=0; dof!=fe.dofs_per_quad; ++dof) + dofs_on_mother[next_index++] = this_face->dof_index(dof, fe_index); + Assert (next_index == dofs_on_mother.size(), + ExcInternalError()); - // dof numbers on the centers of - // the lines bounding this face - for (unsigned int line=0; line<4; ++line) - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - dofs_on_children[next_index++] - = face->line(line)->child(0)->vertex_dof_index(1,dof); + next_index = 0; + + // assert some consistency + // assumptions + Assert ((this_face->child(0)->vertex_index(3) == + this_face->child(1)->vertex_index(2)) && + (this_face->child(0)->vertex_index(3) == + this_face->child(2)->vertex_index(1)) && + (this_face->child(0)->vertex_index(3) == + this_face->child(3)->vertex_index(0)), + ExcInternalError()); + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) + dofs_on_children[next_index++] + = this_face->child(0)->vertex_dof_index(3,dof); - // next the dofs on the lines interior - // to the face; the order of these - // lines is laid down in the - // FiniteElement class documentation - for (unsigned int dof=0; dofchild(0)->line(1)->dof_index(dof); - for (unsigned int dof=0; dofchild(2)->line(1)->dof_index(dof); - for (unsigned int dof=0; dofchild(0)->line(3)->dof_index(dof); - for (unsigned int dof=0; dofchild(1)->line(3)->dof_index(dof); + // dof numbers on the centers of + // the lines bounding this face + for (unsigned int line=0; line<4; ++line) + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) + dofs_on_children[next_index++] + = this_face->line(line)->child(0)->vertex_dof_index(1,dof, fe_index); - // dofs on the bordering lines - for (unsigned int line=0; line<4; ++line) - for (unsigned int child=0; child<2; ++child) - for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) + // next the dofs on the lines interior + // to the face; the order of these + // lines is laid down in the + // FiniteElement class documentation + for (unsigned int dof=0; dofchild(0)->line(1)->dof_index(dof, fe_index); + for (unsigned int dof=0; dofchild(2)->line(1)->dof_index(dof, fe_index); + for (unsigned int dof=0; dofline(line)->child(child)->dof_index(dof); + = this_face->child(0)->line(3)->dof_index(dof, fe_index); + for (unsigned int dof=0; dofchild(1)->line(3)->dof_index(dof, fe_index); - // finally, for the dofs interior - // to the four child faces - for (unsigned int child=0; child<4; ++child) - for (unsigned int dof=0; dof!=fe.dofs_per_quad; ++dof) - dofs_on_children[next_index++] - = face->child(child)->dof_index(dof); - Assert (next_index == dofs_on_children.size(), - ExcInternalError()); + // dofs on the bordering lines + for (unsigned int line=0; line<4; ++line) + for (unsigned int child=0; child<2; ++child) + for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) + dofs_on_children[next_index++] + = this_face->line(line)->child(child)->dof_index(dof, fe_index); - // for each row in the constraint - // matrix for this line: - for (unsigned int row=0; row!=dofs_on_children.size(); ++row) - { - constraints.add_line (dofs_on_children[row]); - for (unsigned int i=0; i!=dofs_on_mother.size(); ++i) - constraints.add_entry (dofs_on_children[row], - dofs_on_mother[i], - fe.constraints()(row,i)); + // finally, for the dofs interior + // to the four child faces + for (unsigned int child=0; child<4; ++child) + for (unsigned int dof=0; dof!=fe.dofs_per_quad; ++dof) + dofs_on_children[next_index++] + = this_face->child(child)->dof_index(dof, fe_index); + Assert (next_index == dofs_on_children.size(), + ExcInternalError()); + + // for each row in the constraint + // matrix for this line: + for (unsigned int row=0; row!=dofs_on_children.size(); ++row) + { + constraints.add_line (dofs_on_children[row]); + for (unsigned int i=0; i!=dofs_on_mother.size(); ++i) + constraints.add_entry (dofs_on_children[row], + dofs_on_mother[i], + fe.constraints()(row,i)); + } } - } + else + { + // this face has no + // children, but it + // could still be that + // it is shared by two + // cells that use a + // different fe index + Assert (cell->face(face)->n_active_fe_indices() == 1, + ExcNotImplemented()); + Assert (cell->face(face) + ->fe_index_is_active(cell->active_fe_index()) == true, + ExcInternalError()); + } + } +#endif + + } } -void DoFTools::make_hanging_node_constraints ( - const hp::DoFHandler<3> &/*dof_handler*/, - ConstraintMatrix &/*constraints*/) + +template +void +DoFTools::make_hanging_node_constraints (const DH &dof_handler, + ConstraintMatrix &constraints) { -//TODO:[?] Implement (required for continuous elements) - Assert (false, ExcNotImplemented()); + // distribute this one template + // function to the appropriate + // functions in 1d, 2d, and 3d + do_make_hanging_node_constraints (dof_handler, + constraints, + internal::DoFTools::int2type()); } -#endif - template @@ -4237,6 +4386,20 @@ DoFTools::make_flux_sparsity_pattern,CompressedBlo #endif +template +void +DoFTools::make_hanging_node_constraints (const DoFHandler &dof_handler, + ConstraintMatrix &constraints); +template +void +DoFTools::make_hanging_node_constraints (const MGDoFHandler &dof_handler, + ConstraintMatrix &constraints); +template +void +DoFTools::make_hanging_node_constraints (const hp::DoFHandler &dof_handler, + ConstraintMatrix &constraints); + + template void -- 2.39.5