From: guido Date: Thu, 6 Jul 2006 14:58:35 +0000 (+0000) Subject: multigrid adapted to level-less faces and edges X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=e746a0ae1a94794042b758ec589116c74705e79c;p=dealii-svn.git multigrid adapted to level-less faces and edges git-svn-id: https://svn.dealii.org/trunk@13348 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/multigrid/mg_dof_accessor.h b/deal.II/deal.II/include/multigrid/mg_dof_accessor.h index 6b52bcc4ac..5c6427ef11 100644 --- a/deal.II/deal.II/include/multigrid/mg_dof_accessor.h +++ b/deal.II/deal.II/include/multigrid/mg_dof_accessor.h @@ -153,7 +153,8 @@ class MGDoFAccessor : public MGDoFObjectAccessor_Inheritance::Ba * vertex for the level this * object lives on. */ - unsigned int mg_vertex_dof_index (const unsigned int vertex, + unsigned int mg_vertex_dof_index (const int level, + const unsigned int vertex, const unsigned int i) const; /** @@ -161,7 +162,8 @@ class MGDoFAccessor : public MGDoFObjectAccessor_Inheritance::Ba * on the @p vertexth vertex to @p index * for the level this object lives on. */ - void set_mg_vertex_dof_index (const unsigned int vertex, + void set_mg_vertex_dof_index (const int level, + const unsigned int vertex, const unsigned int i, const unsigned int index) const; @@ -171,14 +173,16 @@ class MGDoFAccessor : public MGDoFObjectAccessor_Inheritance::Ba * on the level this line lives * on. */ - unsigned int mg_dof_index (const unsigned int i) const; + unsigned int mg_dof_index (const int level, + const unsigned int i) const; /** * Set the index of the @p ith degree * of freedom of this line on the * level this line lives on to @p index. */ - void set_mg_dof_index (const unsigned int i, + void set_mg_dof_index (const int level, + const unsigned int i, const unsigned int index) const; /** @@ -287,7 +291,8 @@ class MGDoFObjectAccessor<1, dim> : public MGDoFAccessor<1,dim> * indices refer to the local numbering * for the level this line lives on. */ - void get_mg_dof_indices (std::vector &dof_indices) const; + void get_mg_dof_indices (const int level, + std::vector &dof_indices) const; /** * Return the value of the given vector @@ -307,7 +312,8 @@ class MGDoFObjectAccessor<1, dim> : public MGDoFAccessor<1,dim> * freedom on this level. */ template - void get_mg_dof_values (const Vector &values, + void get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; }; @@ -370,7 +376,8 @@ class MGDoFObjectAccessor<2, dim> : public MGDoFAccessor<2,dim> * indices refer to the local numbering * for the level this quad lives on. */ - void get_mg_dof_indices (std::vector &dof_indices) const; + void get_mg_dof_indices (const int level, + std::vector &dof_indices) const; /** * Return the value of the given vector @@ -390,7 +397,8 @@ class MGDoFObjectAccessor<2, dim> : public MGDoFAccessor<2,dim> * freedom on this level. */ template - void get_mg_dof_values (const Vector &values, + void get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; /** @@ -461,7 +469,8 @@ class MGDoFObjectAccessor<3, dim> : public MGDoFAccessor<3,dim> * indices refer to the local numbering * for the level this hex lives on. */ - void get_mg_dof_indices (std::vector &dof_indices) const; + void get_mg_dof_indices (const int level, + std::vector &dof_indices) const; /** * Return the value of the given vector @@ -481,7 +490,8 @@ class MGDoFObjectAccessor<3, dim> : public MGDoFAccessor<3,dim> * freedom on this level. */ template - void get_mg_dof_values (const Vector &values, + void get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; /** @@ -540,6 +550,41 @@ class MGDoFCellAccessor : public MGDoFObjectAccessor const int index, const AccessorData *local_data); + /** + * Return the indices of the dofs of this + * hex in the standard ordering: dofs + * on vertex 0, dofs on vertex 1, etc, + * dofs on line 0, dofs on line 1, etc, + * dofs on quad 0, etc. + * + * It is assumed that the vector already + * has the right size beforehand. The + * indices refer to the local numbering + * for the level this hex lives on. + */ + void get_mg_dof_indices (std::vector &dof_indices) const; + + /** + * Return the value of the given vector + * restricted to the dofs of this + * cell in the standard ordering: dofs + * on vertex 0, dofs on vertex 1, etc, + * dofs on line 0, dofs on line 1, etc, + * dofs on quad 0, etc. + * + * It is assumed that the vector already + * has the right size beforehand. The + * indices refer to the multilevel + * numbering local to the present + * level of this cell. The vector shall + * therefore have the same number of + * entries as there are degrees of + * freedom on this level. + */ + template + void get_mg_dof_values (const Vector &values, + Vector &dof_values) const; + /** * Return the @p ith neighbor as a MGDoF cell * iterator. This function is needed since diff --git a/deal.II/deal.II/source/multigrid/mg_dof_accessor.cc b/deal.II/deal.II/source/multigrid/mg_dof_accessor.cc index 943ab81ee9..6500466da9 100644 --- a/deal.II/deal.II/source/multigrid/mg_dof_accessor.cc +++ b/deal.II/deal.II/source/multigrid/mg_dof_accessor.cc @@ -80,7 +80,8 @@ MGDoFAccessor::operator = (const MGDoFAccessor &da) template unsigned int -MGDoFAccessor::mg_vertex_dof_index (const unsigned int vertex, +MGDoFAccessor::mg_vertex_dof_index (const int level, + const unsigned int vertex, const unsigned int i) const { Assert (this->dof_handler != 0, @@ -95,13 +96,14 @@ MGDoFAccessor::mg_vertex_dof_index (const unsigned int vertex, ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_vertex)); return (this->mg_dof_handler->mg_vertex_dofs[this->vertex_index(vertex)] - .get_index (this->present_level, i, this->dof_handler->get_fe().dofs_per_vertex)); + .get_index (level, i, this->dof_handler->get_fe().dofs_per_vertex)); } template void -MGDoFAccessor::set_mg_vertex_dof_index (const unsigned int vertex, +MGDoFAccessor::set_mg_vertex_dof_index (const int level, + const unsigned int vertex, const unsigned int i, const unsigned int index) const { @@ -117,17 +119,18 @@ MGDoFAccessor::set_mg_vertex_dof_index (const unsigned int verte ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_vertex)); this->mg_dof_handler->mg_vertex_dofs[this->vertex_index(vertex)] - .set_index (this->present_level, i, this->dof_handler->get_fe().dofs_per_vertex, index); + .set_index (level, i, this->dof_handler->get_fe().dofs_per_vertex, index); } template unsigned int -MGDoFAccessor::mg_dof_index (const unsigned int i) const +MGDoFAccessor::mg_dof_index (const int level, + const unsigned int i) const { return this->mg_dof_handler - ->template get_dof_index(this->present_level, + ->template get_dof_index(level, this->present_index, 0, i); @@ -135,11 +138,12 @@ MGDoFAccessor::mg_dof_index (const unsigned int i) const template -void MGDoFAccessor::set_mg_dof_index (const unsigned int i, - const unsigned int index) const +void MGDoFAccessor::set_mg_dof_index (const int level, + const unsigned int i, + const unsigned int index) const { this->mg_dof_handler - ->template set_dof_index(this->present_level, + ->template set_dof_index(level, this->present_index, 0, i, @@ -215,7 +219,8 @@ MGDoFObjectAccessor<1, dim>::MGDoFObjectAccessor (const Triangulation *tria template void -MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (std::vector &dof_indices) const +MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (const int level, + std::vector &dof_indices) const { Assert (this->dof_handler != 0, typename BaseClass::ExcInvalidObject()); @@ -232,9 +237,9 @@ MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (std::vector &dof_ std::vector::iterator next = dof_indices.begin(); for (unsigned int vertex=0; vertex<2; ++vertex) for (unsigned int d=0; dmg_vertex_dof_index(vertex,d); + *next++ = this->mg_vertex_dof_index(level,vertex,d); for (unsigned int d=0; dmg_dof_index(d); + *next++ = this->mg_dof_index(level,d); Assert (next == dof_indices.end(), ExcInternalError()); @@ -244,7 +249,8 @@ MGDoFObjectAccessor<1, dim>::get_mg_dof_indices (std::vector &dof_ template template void -MGDoFObjectAccessor<1,dim>::get_mg_dof_values (const Vector &values, +MGDoFObjectAccessor<1,dim>::get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const { Assert (this->dof_handler != 0, @@ -263,9 +269,9 @@ MGDoFObjectAccessor<1,dim>::get_mg_dof_values (const Vector &values, typename Vector::iterator next_dof_value=dof_values.begin(); for (unsigned int vertex=0; vertex<2; ++vertex) for (unsigned int d=0; dmg_vertex_dof_index(vertex,d)); + *next_dof_value++ = values(this->mg_vertex_dof_index(level,vertex,d)); for (unsigned int d=0; dmg_dof_index(d)); + *next_dof_value++ = values(this->mg_dof_index(level,d)); Assert (next_dof_value == dof_values.end(), ExcInternalError()); @@ -289,7 +295,8 @@ MGDoFObjectAccessor<2, dim>::MGDoFObjectAccessor (const Triangulation *tria template void -MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (std::vector &dof_indices) const +MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (const int level, + std::vector &dof_indices) const { Assert (this->dof_handler != 0, typename BaseClass::ExcInvalidObject()); @@ -308,12 +315,12 @@ MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (std::vector &dof_ std::vector::iterator next = dof_indices.begin(); for (unsigned int vertex=0; vertex<4; ++vertex) for (unsigned int d=0; dmg_vertex_dof_index(vertex,d); + *next++ = this->mg_vertex_dof_index(level,vertex,d); for (unsigned int line=0; line<4; ++line) for (unsigned int d=0; dline(line)->mg_dof_index(d); + *next++ = this->line(line)->mg_dof_index(level,d); for (unsigned int d=0; dmg_dof_index(d); + *next++ = this->mg_dof_index(level,d); Assert (next == dof_indices.end(), ExcInternalError()); @@ -323,7 +330,8 @@ MGDoFObjectAccessor<2, dim>::get_mg_dof_indices (std::vector &dof_ template template void -MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const Vector &values, +MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const { Assert (this->dof_handler != 0, @@ -334,7 +342,7 @@ MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const Vector &values, typename BaseClass::ExcInvalidObject()); Assert (dof_values.size() == this->dof_handler->get_fe().dofs_per_cell, typename BaseClass::ExcVectorDoesNotMatch()); - Assert (values.size() == this->mg_dof_handler->n_dofs(this->present_level), + Assert (values.size() == this->mg_dof_handler->n_dofs(level), typename BaseClass::ExcVectorDoesNotMatch()); const unsigned int dofs_per_vertex = this->dof_handler->get_fe().dofs_per_vertex, @@ -343,12 +351,12 @@ MGDoFObjectAccessor<2,dim>::get_mg_dof_values (const Vector &values, typename Vector::iterator next_dof_value=dof_values.begin(); for (unsigned int vertex=0; vertex<4; ++vertex) for (unsigned int d=0; dmg_vertex_dof_index(vertex,d)); + *next_dof_value++ = values(this->mg_vertex_dof_index(level,vertex,d)); for (unsigned int line=0; line<4; ++line) for (unsigned int d=0; dline(line)->mg_dof_index(d)); + *next_dof_value++ = values(this->line(line)->mg_dof_index(level,d)); for (unsigned int d=0; dmg_dof_index(d)); + *next_dof_value++ = values(this->mg_dof_index(level,d)); Assert (next_dof_value == dof_values.end(), ExcInternalError()); @@ -387,7 +395,8 @@ MGDoFObjectAccessor<3, dim>::MGDoFObjectAccessor (const Triangulation *tria template void -MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (std::vector &dof_indices) const +MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (const int level, + std::vector &dof_indices) const { Assert (this->dof_handler != 0, typename BaseClass::ExcInvalidObject()); @@ -408,15 +417,15 @@ MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (std::vector &dof_ std::vector::iterator next = dof_indices.begin(); for (unsigned int vertex=0; vertex<8; ++vertex) for (unsigned int d=0; dmg_vertex_dof_index(vertex,d); + *next++ = this->mg_vertex_dof_index(level,vertex,d); for (unsigned int line=0; line<12; ++line) for (unsigned int d=0; dline(line)->mg_dof_index(d); + *next++ = this->line(line)->mg_dof_index(level,d); for (unsigned int quad=0; quad<6; ++quad) for (unsigned int d=0; dquad(quad)->mg_dof_index(d); + *next++ = this->quad(quad)->mg_dof_index(level,d); for (unsigned int d=0; dmg_dof_index(d); + *next++ = this->mg_dof_index(level,d); Assert (next == dof_indices.end(), ExcInternalError()); @@ -426,7 +435,8 @@ MGDoFObjectAccessor<3, dim>::get_mg_dof_indices (std::vector &dof_ template template void -MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const Vector &values, +MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const { Assert (this->dof_handler != 0, @@ -437,7 +447,7 @@ MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const Vector &values, typename BaseClass::ExcInvalidObject()); Assert (dof_values.size() == this->dof_handler->get_fe().dofs_per_cell, typename BaseClass::ExcVectorDoesNotMatch()); - Assert (values.size() == this->mg_dof_handler->n_dofs(this->present_level), + Assert (values.size() == this->mg_dof_handler->n_dofs(level), typename BaseClass::ExcVectorDoesNotMatch()); const unsigned int dofs_per_vertex = this->dof_handler->get_fe().dofs_per_vertex, @@ -447,13 +457,13 @@ MGDoFObjectAccessor<3,dim>::get_mg_dof_values (const Vector &values, typename Vector::iterator next_dof_value=dof_values.begin(); for (unsigned int vertex=0; vertex<8; ++vertex) for (unsigned int d=0; dmg_vertex_dof_index(vertex,d)); + *next_dof_value++ = values(this->mg_vertex_dof_index(level,vertex,d)); for (unsigned int line=0; line<12; ++line) for (unsigned int d=0; dline(line)->mg_dof_index(d)); + *next_dof_value++ = values(this->line(line)->mg_dof_index(level,d)); for (unsigned int quad=0; quad<6; ++quad) for (unsigned int d=0; dquad(quad)->mg_dof_index(d)); + *next_dof_value++ = values(this->quad(quad)->mg_dof_index(level,d)); for (unsigned int d=0; ddof_index(d)); @@ -507,6 +517,25 @@ MGDoFCellAccessor::MGDoFCellAccessor (const Triangulation *tria, {} +template +void +MGDoFCellAccessor::get_mg_dof_indices (std::vector &dof_indices) const +{ + MGDoFObjectAccessor::get_mg_dof_indices (this->level(), dof_indices); +} + + +template +template +void +MGDoFCellAccessor::get_mg_dof_values (const Vector &values, + Vector &dof_values) const +{ + MGDoFObjectAccessor::get_mg_dof_values (this->level(), values, dof_values); +} + + + template TriaIterator > MGDoFCellAccessor::neighbor (const unsigned int i) const @@ -601,13 +630,15 @@ neighbor_child_on_subface (const unsigned int face, template void MGDoFObjectAccessor<1,deal_II_dimension>:: -get_mg_dof_values (const Vector &values, +get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; template void MGDoFObjectAccessor<1,deal_II_dimension>:: -get_mg_dof_values (const Vector &values, +get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; @@ -616,13 +647,15 @@ get_mg_dof_values (const Vector &values, template void MGDoFObjectAccessor<2,deal_II_dimension>:: -get_mg_dof_values (const Vector &values, +get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; template void MGDoFObjectAccessor<2,deal_II_dimension>:: -get_mg_dof_values (const Vector &values, +get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; #endif @@ -633,16 +666,30 @@ get_mg_dof_values (const Vector &values, template void MGDoFObjectAccessor<3,deal_II_dimension>:: -get_mg_dof_values (const Vector &values, +get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; template void MGDoFObjectAccessor<3,deal_II_dimension>:: -get_mg_dof_values (const Vector &values, +get_mg_dof_values (const int level, + const Vector &values, Vector &dof_values) const; #endif +template +void +MGDoFCellAccessor:: +get_mg_dof_values (const Vector &values, + Vector &dof_values) const; + +template +void +MGDoFCellAccessor:: +get_mg_dof_values (const Vector &values, + Vector &dof_values) const; + #if deal_II_dimension == 1 diff --git a/deal.II/deal.II/source/multigrid/mg_dof_handler.cc b/deal.II/deal.II/source/multigrid/mg_dof_handler.cc index bbd5a8b1c2..d4c699166c 100644 --- a/deal.II/deal.II/source/multigrid/mg_dof_handler.cc +++ b/deal.II/deal.II/source/multigrid/mg_dof_handler.cc @@ -1733,12 +1733,12 @@ MGDoFHandler<1>::distribute_dofs_on_cell (cell_iterator &cell, { if (v==0) for (unsigned int d=0; dselected_fe->dofs_per_vertex; ++d) - cell->set_mg_vertex_dof_index (0, d, - neighbor->mg_vertex_dof_index (1, d)); + cell->set_mg_vertex_dof_index (cell->level(), 0, d, + neighbor->mg_vertex_dof_index (cell->level(), 1, d)); else for (unsigned int d=0; dselected_fe->dofs_per_vertex; ++d) - cell->set_mg_vertex_dof_index (1, d, - neighbor->mg_vertex_dof_index (0, d)); + cell->set_mg_vertex_dof_index (cell->level(), 1, d, + neighbor->mg_vertex_dof_index (cell->level(), 0, d)); // next neighbor continue; @@ -1747,13 +1747,13 @@ MGDoFHandler<1>::distribute_dofs_on_cell (cell_iterator &cell, // otherwise: create dofs newly for (unsigned int d=0; dselected_fe->dofs_per_vertex; ++d) - cell->set_mg_vertex_dof_index (v, d, next_free_dof++); + cell->set_mg_vertex_dof_index (cell->level(), v, d, next_free_dof++); }; // dofs of line if (this->selected_fe->dofs_per_line > 0) for (unsigned int d=0; dselected_fe->dofs_per_line; ++d) - cell->set_mg_dof_index (d, next_free_dof++); + cell->set_mg_dof_index (cell->level(), d, next_free_dof++); // note that this cell has been processed cell->set_user_flag (); @@ -1776,9 +1776,9 @@ MGDoFHandler<2>::distribute_dofs_on_cell (cell_iterator &cell, // check whether dofs for this // vertex have been distributed // (only check the first dof) - if (cell->mg_vertex_dof_index(vertex, 0) == DoFHandler<2>::invalid_dof_index) + if (cell->mg_vertex_dof_index(cell->level(), vertex, 0) == DoFHandler<2>::invalid_dof_index) for (unsigned int d=0; dselected_fe->dofs_per_vertex; ++d) - cell->set_mg_vertex_dof_index (vertex, d, next_free_dof++); + cell->set_mg_vertex_dof_index (cell->level(), vertex, d, next_free_dof++); // for the four sides if (this->selected_fe->dofs_per_line > 0) @@ -1789,17 +1789,17 @@ MGDoFHandler<2>::distribute_dofs_on_cell (cell_iterator &cell, // distribute dofs if necessary: // check whether line dof is already // numbered (check only first dof) - if (line->mg_dof_index(0) == DoFHandler<2>::invalid_dof_index) + if (line->mg_dof_index(cell->level(), 0) == DoFHandler<2>::invalid_dof_index) // if not: distribute dofs for (unsigned int d=0; dselected_fe->dofs_per_line; ++d) - line->set_mg_dof_index (d, next_free_dof++); + line->set_mg_dof_index (cell->level(), d, next_free_dof++); }; // dofs of quad if (this->selected_fe->dofs_per_quad > 0) for (unsigned int d=0; dselected_fe->dofs_per_quad; ++d) - cell->set_mg_dof_index (d, next_free_dof++); + cell->set_mg_dof_index (cell->level(), d, next_free_dof++); // note that this cell has been processed @@ -1823,9 +1823,9 @@ MGDoFHandler<3>::distribute_dofs_on_cell (cell_iterator &cell, // check whether dofs for this // vertex have been distributed // (only check the first dof) - if (cell->mg_vertex_dof_index(vertex, 0) == DoFHandler<3>::invalid_dof_index) + if (cell->mg_vertex_dof_index(cell->level(), vertex, 0) == DoFHandler<3>::invalid_dof_index) for (unsigned int d=0; dselected_fe->dofs_per_vertex; ++d) - cell->set_mg_vertex_dof_index (vertex, d, next_free_dof++); + cell->set_mg_vertex_dof_index (cell->level(), vertex, d, next_free_dof++); // for the lines if (this->selected_fe->dofs_per_line > 0) @@ -1836,10 +1836,10 @@ MGDoFHandler<3>::distribute_dofs_on_cell (cell_iterator &cell, // distribute dofs if necessary: // check whether line dof is already // numbered (check only first dof) - if (line->mg_dof_index(0) == DoFHandler<3>::invalid_dof_index) + if (line->mg_dof_index(cell->level(), 0) == DoFHandler<3>::invalid_dof_index) // if not: distribute dofs for (unsigned int d=0; dselected_fe->dofs_per_line; ++d) - line->set_mg_dof_index (d, next_free_dof++); + line->set_mg_dof_index (cell->level(), d, next_free_dof++); }; // for the quads @@ -1851,17 +1851,17 @@ MGDoFHandler<3>::distribute_dofs_on_cell (cell_iterator &cell, // distribute dofs if necessary: // check whether line dof is already // numbered (check only first dof) - if (quad->mg_dof_index(0) == DoFHandler<3>::invalid_dof_index) + if (quad->mg_dof_index(cell->level(), 0) == DoFHandler<3>::invalid_dof_index) // if not: distribute dofs for (unsigned int d=0; dselected_fe->dofs_per_quad; ++d) - quad->set_mg_dof_index (d, next_free_dof++); + quad->set_mg_dof_index (cell->level(), d, next_free_dof++); }; // dofs of cell if (this->selected_fe->dofs_per_hex > 0) for (unsigned int d=0; dselected_fe->dofs_per_hex; ++d) - cell->set_mg_dof_index (d, next_free_dof++); + cell->set_mg_dof_index (cell->level(), d, next_free_dof++); // note that this cell has been processed @@ -1975,7 +1975,7 @@ void MGDoFHandler<2>::renumber_dofs (const unsigned int level, if (line->user_flag_set()) { for (unsigned int d=0; dselected_fe->dofs_per_line; ++d) - line->set_mg_dof_index (d, new_numbers[line->mg_dof_index(d)]); + line->set_mg_dof_index (level, d, new_numbers[line->mg_dof_index(level, d)]); line->clear_user_flag(); } // finally, restore user flags @@ -2042,7 +2042,7 @@ void MGDoFHandler<3>::renumber_dofs (const unsigned int level, if (line->user_flag_set()) { for (unsigned int d=0; dselected_fe->dofs_per_line; ++d) - line->set_mg_dof_index (d, new_numbers[line->mg_dof_index(d)]); + line->set_mg_dof_index (level, d, new_numbers[line->mg_dof_index(level, d)]); line->clear_user_flag(); } // finally, restore user flags @@ -2074,7 +2074,7 @@ void MGDoFHandler<3>::renumber_dofs (const unsigned int level, if (quad->user_flag_set()) { for (unsigned int d=0; dselected_fe->dofs_per_quad; ++d) - quad->set_mg_dof_index (d, new_numbers[quad->mg_dof_index(d)]); + quad->set_mg_dof_index (level, d, new_numbers[quad->mg_dof_index(level, d)]); quad->clear_user_flag(); } // finally, restore user flags diff --git a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc index 124e7b3450..b10aca91d4 100644 --- a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc +++ b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc @@ -1111,7 +1111,7 @@ MGTools::make_boundary_list( // boundary values of dofs on this // face face_dofs.resize (fe.dofs_per_face); - face->get_mg_dof_indices (face_dofs); + face->get_mg_dof_indices (level, face_dofs); if (fe_is_system) { // enter those dofs diff --git a/deal.II/deal.II/source/multigrid/mg_smoother.cc b/deal.II/deal.II/source/multigrid/mg_smoother.cc index f71d9e9875..e916067704 100644 --- a/deal.II/deal.II/source/multigrid/mg_smoother.cc +++ b/deal.II/deal.II/source/multigrid/mg_smoother.cc @@ -91,7 +91,7 @@ MGSmootherContinuous::MGSmootherContinuous ( == level-1)) { // get indices of this face - cell->face(face)->get_mg_dof_indices (dofs_on_face); + cell->face(face)->get_mg_dof_indices (level, dofs_on_face); // append them to the levelwise // list boundary_dofs.insert (boundary_dofs.end(),