From: Wolfgang Bangerth Date: Mon, 15 May 2000 12:04:17 +0000 (+0000) Subject: Don't use the user flags any more for make_hanging_node_constraints. X-Git-Tag: v8.0.0~20550 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7eea11f59384da2e6d5fb45621a3fb972a8a300e;p=dealii.git Don't use the user flags any more for make_hanging_node_constraints. git-svn-id: https://svn.dealii.org/trunk@2855 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 80b2b24b2b..785142f226 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -318,6 +318,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler &dof, } + #if deal_II_dimension == 1 template <> @@ -330,6 +331,7 @@ void DoFTools::make_hanging_node_constraints (const DoFHandler<1> &, #endif + #if deal_II_dimension == 2 template <> @@ -337,91 +339,87 @@ void DoFTools::make_hanging_node_constraints (const DoFHandler<2> &dof_handler, ConstraintMatrix &constraints) { const unsigned int dim = 2; - const Triangulation &tria = dof_handler.get_tria(); - const FiniteElement &fe = dof_handler.get_fe(); - // first mark all faces which are subject - // to 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 - DoFHandler::line_iterator line = dof_handler.begin_line(), - endl = dof_handler.end_line(); - for (; line!=endl; ++line) - line->clear_user_flag (); + const FiniteElement &fe = dof_handler.get_fe(); - Triangulation::active_cell_iterator cell = tria.begin_active(), - endc = tria.end(); - for (; cell!=endc; ++cell) - for (unsigned int face=0; face::faces_per_cell; ++face) - if (cell->face(face)->has_children()) - cell->face(face)->set_user_flag(); - // have space for the degrees of // freedom on mother and child // lines - vector dofs_on_mother(2*fe.dofs_per_vertex+ - fe.dofs_per_line); - vector dofs_on_children(fe.dofs_per_vertex+ - 2*fe.dofs_per_line); - - Assert(2*fe.dofs_per_vertex+fe.dofs_per_line == - fe.constraints().n(), - ExcDimensionMismatch(2*fe.dofs_per_vertex+ - fe.dofs_per_line, + 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; + + vector dofs_on_mother(n_dofs_on_mother); + 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(fe.dofs_per_vertex+2*fe.dofs_per_line == - fe.constraints().m(), - ExcDimensionMismatch(3*fe.dofs_per_vertex+ - 2*fe.dofs_per_line, + 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. - for (line = dof_handler.begin_line(); line != endl; ++line) - // if dofs on this line are subject - // to constraints - if (line->user_flag_set() == true) - { - // fill the dofs indices. Use same - // enumeration scheme as in - // #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) + // 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 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 + // #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_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) - { - 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)); - }; - }; + 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) + { + 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)); + }; + }; }; #endif + #if deal_II_dimension == 3 template <> @@ -429,152 +427,138 @@ void DoFTools::make_hanging_node_constraints (const DoFHandler<3> &dof_handler, ConstraintMatrix &constraints) { const unsigned int dim = 3; - const Triangulation &tria = dof_handler.get_tria(); + const FiniteElement &fe = dof_handler.get_fe(); - // first mark all faces which are subject - // to 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 - DoFHandler::face_iterator face = dof_handler.begin_face(), - endf = dof_handler.end_face(); - for (; face!=endf; ++face) - face->clear_user_flag (); + // 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); + + vector dofs_on_mother(n_dofs_on_mother); + 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())); - Triangulation::active_cell_iterator cell = tria.begin_active(), - endc = tria.end(); + // 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 face=0; face::faces_per_cell; ++face) if (cell->face(face)->has_children()) - cell->face(face)->set_user_flag(); - - - // allocate space for dof indices - vector dofs_on_mother(4*fe.dofs_per_vertex+ - 4*fe.dofs_per_line+ - fe.dofs_per_quad); - vector dofs_on_children(5*fe.dofs_per_vertex+ - 12*fe.dofs_per_line+ - 4*fe.dofs_per_quad); - Assert(4*fe.dofs_per_vertex+ - 4*fe.dofs_per_line+ - fe.dofs_per_quad - == - fe.constraints().n(), - ExcDimensionMismatch(4*fe.dofs_per_vertex+ - 4*fe.dofs_per_line+ - fe.dofs_per_quad, - fe.constraints().n())); - Assert(5*fe.dofs_per_vertex+ - 12*fe.dofs_per_line+ - 4*fe.dofs_per_quad - == - fe.constraints().m(), - ExcDimensionMismatch(5*fe.dofs_per_vertex+ - 12*fe.dofs_per_line+ - 4*fe.dofs_per_quad, - fe.constraints().m())); - - // loop over all faces; only on faces - // there can be constraints. - for (face=dof_handler.begin_face(); face != endf; ++face) - // if dofs on this line are subject - // to constraints - if (face->user_flag_set() == true) - { - 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()); - - // dof numbers on vertex at the center - // of the face, which is vertex 2 of - // child zero, or vertex 3 of child 1 - // or vertex 0 of child 2 or vertex 1 - // of child 3. We're a bit cautious and - // check this (also an additional safety - // check for the internal states of the - // library) - Assert ((face->child(0)->vertex_dof_index(2,0) == - face->child(1)->vertex_dof_index(3,0)) && - (face->child(0)->vertex_dof_index(2,0) == - face->child(2)->vertex_dof_index(0,0)) && - (face->child(0)->vertex_dof_index(2,0) == - face->child(3)->vertex_dof_index(1,0)), - ExcInternalError()); - next_index = 0; - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - dofs_on_children[next_index++] = face->child(0)->vertex_dof_index(2,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++] = face->line(line)->child(0)->vertex_dof_index(1,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(1)->line(2)->dof_index(dof); - for (unsigned int dof=0; dofchild(2)->line(3)->dof_index(dof); - for (unsigned int dof=0; dofchild(3)->line(0)->dof_index(dof); - - // dofs on the bordering lines - for (unsigned int line=0; line<4; ++line) - for (unsigned int child=0; child<2; ++child) + { + const DoFHandler::face_iterator face = cell->face(face); + + // fill the dofs indices. Use same + // enumeration scheme as in + // #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_children[next_index++] = face->line(line)->child(child)->dof_index(dof); - - // finally, for the dofs interior - // to the four child faces - for (unsigned int child=0; child<4; ++child) + dofs_on_mother[next_index++] = face->line(line)->dof_index(dof); 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()); - - Assert (dofs_on_children.size() == - fe.constraints().m(), - ExcDimensionMismatch(dofs_on_children.size(), - fe.constraints().m())); - Assert (dofs_on_mother.size() == - fe.constraints().n(), - ExcDimensionMismatch(dofs_on_mother.size(), - fe.constraints().n())); - - // 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)); - }; - }; + dofs_on_mother[next_index++] = face->dof_index(dof); + Assert (next_index == dofs_on_mother.size(), + ExcInternalError()); + + next_index = 0; + Assert ((face->child(0)->vertex_dof_index(2,0) == + face->child(1)->vertex_dof_index(3,0)) && + (face->child(0)->vertex_dof_index(2,0) == + face->child(2)->vertex_dof_index(0,0)) && + (face->child(0)->vertex_dof_index(2,0) == + face->child(3)->vertex_dof_index(1,0)), + ExcInternalError()); + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) + dofs_on_children[next_index++] + = face->child(0)->vertex_dof_index(2,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++] + = face->line(line)->child(0)->vertex_dof_index(1,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(1)->line(2)->dof_index(dof); + for (unsigned int dof=0; dofchild(2)->line(3)->dof_index(dof); + for (unsigned int dof=0; dofchild(3)->line(0)->dof_index(dof); + + // 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++] + = face->line(line)->child(child)->dof_index(dof); + + // 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()); + + // 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)); + }; + }; }; #endif + template void DoFTools::distribute_cell_to_dof_vector (const DoFHandler &dof_handler, const Vector &cell_data, @@ -818,42 +802,26 @@ DoFTools::extract_hanging_node_dofs (const DoFHandler<2> &dof_handler, // preset all values by false fill_n (selected_dofs.begin(), dof_handler.n_dofs(), false); - const Triangulation &tria = dof_handler.get_tria(); const FiniteElement &fe = dof_handler.get_fe(); - - // first mark all faces which are subject - // to 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 - DoFHandler::line_iterator line = dof_handler.begin_line(), - endl = dof_handler.end_line(); - for (; line!=endl; ++line) - line->clear_user_flag (); - - Triangulation::active_cell_iterator cell = tria.begin_active(), - endc = tria.end(); + + // this function is similar to the + // make_sparsity_pattern function, + // see there for more information + DoFHandler::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()) - cell->face(face)->set_user_flag(); - - // loop over all lines; only on lines - // can there be constraints. - for (line = dof_handler.begin_line(); line != endl; ++line) - // if dofs on this line are subject - // to constraints - if (line->user_flag_set() == true) - { - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - selected_dofs[line->child(0)->vertex_dof_index(1,dof)] = true; - - for (unsigned int child=0; child<2; ++child) - for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) - selected_dofs[line->child(child)->dof_index(dof)] = true; - }; + { + const DoFHandler::line_iterator line = cell->face(face); + + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) + selected_dofs[line->child(0)->vertex_dof_index(1,dof)] = true; + + for (unsigned int child=0; child<2; ++child) + for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) + selected_dofs[line->child(child)->dof_index(dof)] = true; + }; }; #endif @@ -868,69 +836,60 @@ DoFTools::extract_hanging_node_dofs (const DoFHandler<3> &dof_handler, vector &selected_dofs) { const unsigned int dim = 3; - const Triangulation &tria = dof_handler.get_tria(); + + Assert(selected_dofs.size() == dof_handler.n_dofs(), + ExcDimensionMismatch(selected_dofs.size(), dof_handler.n_dofs())); + // preset all values by false + fill_n (selected_dofs.begin(), dof_handler.n_dofs(), false); + const FiniteElement &fe = dof_handler.get_fe(); - // first mark all faces which are subject - // to 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 - DoFHandler::face_iterator face = dof_handler.begin_face(), - endf = dof_handler.end_face(); - for (; face!=endf; ++face) - face->clear_user_flag (); + // this function is similar to the + // make_sparsity_pattern function, + // see there for more information - Triangulation::active_cell_iterator cell = tria.begin_active(), - endc = tria.end(); + DoFHandler::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()) - cell->face(face)->set_user_flag(); - - // loop over all faces; only on faces - // can there be constraints. - for (face=dof_handler.begin_face(); face != endf; ++face) - // if dofs on this line are subject - // to constraints - if (face->user_flag_set() == true) - { - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - selected_dofs[face->child(0)->vertex_dof_index(2,dof)] = true; - - // dof numbers on the centers of - // the lines bounding this face - for (unsigned int line=0; line<4; ++line) + { + const DoFHandler::face_iterator face = cell->face(face); + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - selected_dofs[face->line(line)->child(0)->vertex_dof_index(1,dof)] = true; - - // 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)] = true; - for (unsigned int dof=0; dofchild(1)->line(2)->dof_index(dof)] = true; - for (unsigned int dof=0; dofchild(2)->line(3)->dof_index(dof)] = true; - for (unsigned int dof=0; dofchild(3)->line(0)->dof_index(dof)] = true; - - // 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) - selected_dofs[face->line(line)->child(child)->dof_index(dof)] = true; + selected_dofs[face->child(0)->vertex_dof_index(2,dof)] = true; + + // 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) + selected_dofs[face->line(line)->child(0)->vertex_dof_index(1,dof)] = true; + + // 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)] = true; + for (unsigned int dof=0; dofchild(1)->line(2)->dof_index(dof)] = true; + for (unsigned int dof=0; dofchild(2)->line(3)->dof_index(dof)] = true; + for (unsigned int dof=0; dofchild(3)->line(0)->dof_index(dof)] = true; + + // 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) + selected_dofs[face->line(line)->child(child)->dof_index(dof)] = true; - // 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) - selected_dofs[face->child(child)->dof_index(dof)] = true; - }; + // 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) + selected_dofs[face->child(child)->dof_index(dof)] = true; + }; }; #endif