From e0c362a15d287399b9f310a2d095e2560a317c80 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Fri, 28 Apr 2006 20:20:49 +0000 Subject: [PATCH] A few more functions to access hp dofs git-svn-id: https://svn.dealii.org/trunk@12938 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/dofs/dof_accessor.templates.h | 458 ++++------------ deal.II/deal.II/include/dofs/hp_dof_levels.h | 63 ++- deal.II/deal.II/source/dofs/dof_accessor.cc | 102 ++-- deal.II/deal.II/source/dofs/hp_dof_levels.cc | 489 +++++++++++++++++- 4 files changed, 675 insertions(+), 437 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_accessor.templates.h b/deal.II/deal.II/include/dofs/dof_accessor.templates.h index 6bcd68b5fc..c31e55eb65 100644 --- a/deal.II/deal.II/include/dofs/dof_accessor.templates.h +++ b/deal.II/deal.II/include/dofs/dof_accessor.templates.h @@ -218,7 +218,7 @@ DoFObjectAccessor<2,DH>::line (const unsigned int i) const this->present_level, this->line_index (i), this->dof_handler - ); + ); } @@ -309,7 +309,7 @@ DoFObjectAccessor<3,DH>::line (const unsigned int i) const this->present_level, l->index(), this->dof_handler - ); + ); } @@ -326,7 +326,7 @@ DoFObjectAccessor<3,DH>::quad (const unsigned int i) const this->present_level, this->quad_index (i), this->dof_handler - ); + ); } @@ -966,10 +966,10 @@ DoFObjectAccessor<1,hp::DoFHandler<1> >::vertex_dof_index (const unsigned int ve const unsigned int fe_index) const { return internal::hp::DoFLevel<0>:: - get_hp_vertex_dof_index (*this->dof_handler, - this->vertex_index(vertex), - fe_index, - i); + get_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i); } @@ -982,10 +982,10 @@ DoFObjectAccessor<1,hp::DoFHandler<2> >::vertex_dof_index (const unsigned int ve const unsigned int fe_index) const { return internal::hp::DoFLevel<0>:: - get_hp_vertex_dof_index (*this->dof_handler, - this->vertex_index(vertex), - fe_index, - i); + get_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i); } @@ -997,10 +997,10 @@ DoFObjectAccessor<1,hp::DoFHandler<3> >::vertex_dof_index (const unsigned int ve const unsigned int fe_index) const { return internal::hp::DoFLevel<0>:: - get_hp_vertex_dof_index (*this->dof_handler, - this->vertex_index(vertex), - fe_index, - i); + get_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i); } @@ -1012,10 +1012,10 @@ DoFObjectAccessor<2,hp::DoFHandler<2> >::vertex_dof_index (const unsigned int ve const unsigned int fe_index) const { return internal::hp::DoFLevel<0>:: - get_hp_vertex_dof_index (*this->dof_handler, - this->vertex_index(vertex), - fe_index, - i); + get_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i); } @@ -1028,10 +1028,10 @@ DoFObjectAccessor<2,hp::DoFHandler<3> >::vertex_dof_index (const unsigned int ve const unsigned int fe_index) const { return internal::hp::DoFLevel<0>:: - get_hp_vertex_dof_index (*this->dof_handler, - this->vertex_index(vertex), - fe_index, - i); + get_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i); } @@ -1044,10 +1044,10 @@ DoFObjectAccessor<3,hp::DoFHandler<3> >::vertex_dof_index (const unsigned int ve const unsigned int fe_index) const { return internal::hp::DoFLevel<0>:: - get_hp_vertex_dof_index (*this->dof_handler, - this->vertex_index(vertex), - fe_index, - i); + get_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i); } @@ -1059,42 +1059,11 @@ unsigned int DoFObjectAccessor<1,hp::DoFHandler<1> >::dof_index (const unsigned int i, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<1>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - - // since the exception classes are - // from a template dependent base - // class, we have to fully qualify - // them. to work around more - // trouble, typedef the template - // dependent base class to a - // non-template dependent name and - // use that to specify the - // qualified exception names - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_line, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_line)); - - // for lines in 1d, we only have to - // take care of the offset, no - // linked lists or whatever (see - // hp::DoFLevels) - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_line_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->line_dofs.size(), - ExcInternalError()); - - return this->dof_handler->levels[this->present_level]->line_dofs[offset+i]; + return this->dof_handler->levels[this->present_level] + ->get_line_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i); } @@ -1144,37 +1113,11 @@ unsigned int DoFObjectAccessor<1,hp::DoFHandler<2> >::dof_index (const unsigned int i, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<2>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - // since the exception classes are - // from a template dependent base - // class, we have to fully qualify - // them. to work around more - // trouble, typedef the template - // dependent base class to a - // non-template dependent name and - // use that to specify the - // qualified exception names - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_line, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_line)); - - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_line_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->line_dofs.size(), - ExcInternalError()); - - return this->dof_handler->levels[this->present_level]->line_dofs[offset+i]; + return this->dof_handler->levels[this->present_level] + ->get_line_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i); } @@ -1224,37 +1167,11 @@ unsigned int DoFObjectAccessor<1,hp::DoFHandler<3> >::dof_index (const unsigned int i, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - // since the exception classes are - // from a template dependent base - // class, we have to fully qualify - // them. to work around more - // trouble, typedef the template - // dependent base class to a - // non-template dependent name and - // use that to specify the - // qualified exception names - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_line, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_line)); - - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_line_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->line_dofs.size(), - ExcInternalError()); - - return this->dof_handler->levels[this->present_level]->line_dofs[offset+i]; + return this->dof_handler->levels[this->present_level] + ->get_line_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i); } @@ -1305,34 +1222,12 @@ DoFObjectAccessor<1,hp::DoFHandler<1> >::set_dof_index (const unsigned int i, const unsigned int index, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<1>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_line, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_line)); - - // for lines in 1d, we only have to - // take care of the offset, no - // linked lists or whatever (see - // hp::DoFLevels) - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_line_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->line_dofs.size(), - ExcInternalError()); - - this->dof_handler->levels[this->present_level]->line_dofs[offset+i] = index; + this->dof_handler->levels[this->present_level] + ->set_line_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i, + index); } @@ -1344,34 +1239,12 @@ DoFObjectAccessor<1, hp::DoFHandler<2> >::set_dof_index (const unsigned int i, const unsigned int index, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<2>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_line, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_line)); - -//TODO:[?] In two dimension it could happen that we have different active_fe_indices -// on a line between to cells. Hence we have to differentiate between these two cases. -// Unfortunately, this requires more information then available now. - - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_line_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->line_dofs.size(), - ExcInternalError()); - - this->dof_handler->levels[this->present_level]->line_dofs[offset+i] = index; + this->dof_handler->levels[this->present_level] + ->set_line_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i, + index); } @@ -1383,30 +1256,12 @@ DoFObjectAccessor<1, hp::DoFHandler<3> >::set_dof_index (const unsigned int i, const unsigned int index, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_line, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_line)); - - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_line_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->line_dofs.size(), - ExcInternalError()); - - this->dof_handler->levels[this->present_level]->line_dofs[offset+i] = index; + this->dof_handler->levels[this->present_level] + ->set_line_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i, + index); } @@ -1419,34 +1274,11 @@ inline unsigned int DoFObjectAccessor<2,hp::DoFHandler<2> >::dof_index (const unsigned int i, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<2>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_quad, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_quad)); - - // for quads in 2d, we only have to - // take care of the offset, no - // linked lists or whatever (see - // hp::DoFLevels) - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_quad_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->quad_dofs.size(), - ExcInternalError()); - - return this->dof_handler->levels[this->present_level]->quad_dofs[offset+i]; + return this->dof_handler->levels[this->present_level] + ->get_quad_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i); } @@ -1456,30 +1288,11 @@ inline unsigned int DoFObjectAccessor<2,hp::DoFHandler<3> >::dof_index (const unsigned int i, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_quad, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_quad)); - - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_quad_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->quad_dofs.size(), - ExcInternalError()); - - return this->dof_handler->levels[this->present_level]->quad_dofs[offset+i]; + return this->dof_handler->levels[this->present_level] + ->get_quad_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i); } @@ -1491,34 +1304,12 @@ DoFObjectAccessor<2, hp::DoFHandler<2> >::set_dof_index (const unsigned int i, const unsigned int index, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<2>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_quad, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_quad)); - - // for quads in 2d, we only have to - // take care of the offset, no - // linked lists or whatever (see - // hp::DoFLevels) - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_quad_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->quad_dofs.size(), - ExcInternalError()); - - this->dof_handler->levels[this->present_level]->quad_dofs[offset+i] = index; + this->dof_handler->levels[this->present_level] + ->set_quad_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i, + index); } @@ -1530,30 +1321,12 @@ DoFObjectAccessor<2, hp::DoFHandler<3> >::set_dof_index (const unsigned int i, const unsigned int index, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_quad, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_quad)); - - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_quad_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->quad_dofs.size(), - ExcInternalError()); - - this->dof_handler->levels[this->present_level]->quad_dofs[offset+i] = index; + this->dof_handler->levels[this->present_level] + ->set_quad_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i, + index); } @@ -1660,34 +1433,11 @@ unsigned int DoFObjectAccessor<3,hp::DoFHandler<3> >::dof_index (const unsigned int i, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_hex, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_hex)); - - // for hexes in 3d, we only have to - // take care of the offset, no - // linked lists or whatever (see - // hp::DoFLevels) - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_hex_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->hex_dofs.size(), - ExcInternalError()); - - return this->dof_handler->levels[this->present_level]->hex_dofs[offset+i]; + return this->dof_handler->levels[this->present_level] + ->get_hex_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i); } @@ -1699,34 +1449,12 @@ DoFObjectAccessor<3, hp::DoFHandler<3> >::set_dof_index (const unsigned int i, const unsigned int index, const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - typedef DoFAccessor > BaseClass; - - Assert (this->dof_handler != 0, - BaseClass::ExcInvalidObject()); - Assert (static_cast(this->present_level) < this->dof_handler->levels.size(), - ExcMessage ("DoFHandler not initialized")); - // make sure a FE has been selected - // and enough room was reserved - Assert (&this->dof_handler->get_fe() != 0, - BaseClass::ExcInvalidObject()); - Assert (idof_handler->get_fe()[fe_index].dofs_per_hex, - ExcIndexRange (i, 0, this->dof_handler->get_fe()[fe_index].dofs_per_hex)); - - // for hexes in 3d, we only have to - // take care of the offset, no - // linked lists or whatever (see - // hp::DoFLevels) - const unsigned int offset = this->dof_handler->levels[this->present_level] - ->dof_hex_index_offset[this->present_index]; - Assert (offset != this->dof_handler->invalid_dof_index, - ExcInternalError()); - Assert (offset+i < - this->dof_handler->levels[this->present_level]->hex_dofs.size(), - ExcInternalError()); - - this->dof_handler->levels[this->present_level]->hex_dofs[offset+i] = index; + this->dof_handler->levels[this->present_level] + ->set_hex_dof_index (*this->dof_handler, + this->present_index, + fe_index, + i, + index); } diff --git a/deal.II/deal.II/include/dofs/hp_dof_levels.h b/deal.II/deal.II/include/dofs/hp_dof_levels.h index db47306c48..3b8f6b8dea 100644 --- a/deal.II/deal.II/include/dofs/hp_dof_levels.h +++ b/deal.II/deal.II/include/dofs/hp_dof_levels.h @@ -89,19 +89,19 @@ namespace internal template static void - set_hp_vertex_dof_index (::hp::DoFHandler &dof_handler, - const unsigned int vertex_index, - const unsigned int fe_index, - const unsigned int local_index, - const unsigned int global_index); + set_vertex_dof_index (::hp::DoFHandler &dof_handler, + const unsigned int vertex_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index); template static unsigned int - get_hp_vertex_dof_index (const ::hp::DoFHandler &dof_handler, - const unsigned int vertex_index, - const unsigned int fe_index, - const unsigned int local_index); + get_vertex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int vertex_index, + const unsigned int fe_index, + const unsigned int local_index); }; @@ -179,6 +179,21 @@ namespace internal */ std::vector line_dofs; + template + void + set_line_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int line_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index); + + template + unsigned int + get_line_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int line_index, + const unsigned int fe_index, + const unsigned int local_index) const; + /** * Determine an estimate for the * memory consumption (in bytes) @@ -216,6 +231,21 @@ namespace internal */ std::vector quad_dofs; + template + void + set_quad_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int quad_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index); + + template + unsigned int + get_quad_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int quad_index, + const unsigned int fe_index, + const unsigned int local_index) const; + /** * Determine an estimate for the * memory consumption (in bytes) @@ -254,6 +284,21 @@ namespace internal */ std::vector hex_dofs; + template + void + set_hex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int hex_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index); + + template + unsigned int + get_hex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int hex_index, + const unsigned int fe_index, + const unsigned int local_index) const; + /** * Determine an estimate for the * memory consumption (in bytes) diff --git a/deal.II/deal.II/source/dofs/dof_accessor.cc b/deal.II/deal.II/source/dofs/dof_accessor.cc index 68528722e5..324e121d3a 100644 --- a/deal.II/deal.II/source/dofs/dof_accessor.cc +++ b/deal.II/deal.II/source/dofs/dof_accessor.cc @@ -65,7 +65,7 @@ template void DoFObjectAccessor<1, DH>::set_vertex_dof_index (const unsigned int vertex, const unsigned int i, const unsigned int index, - const unsigned int fe_index) const + const unsigned int fe_index) const { Assert (fe_index == DoFHandler<1>::default_fe_index, @@ -199,7 +199,7 @@ void DoFObjectAccessor<2, DH>::set_vertex_dof_index (const unsigned int vertex, const unsigned int i, const unsigned int index, - const unsigned int fe_index) const + const unsigned int fe_index) const { Assert (fe_index == DoFHandler<1>::default_fe_index, ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects")); @@ -331,7 +331,7 @@ template void DoFObjectAccessor<3, DH>::set_vertex_dof_index (const unsigned int vertex, const unsigned int i, const unsigned int index, - const unsigned int fe_index) const + const unsigned int fe_index) const { Assert (fe_index == DoFHandler<1>::default_fe_index, ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects")); @@ -448,14 +448,16 @@ DoFObjectAccessor<3,DH>::set_dof_values (const Vector &local_values, #if deal_II_dimension == 1 template <> -void DoFObjectAccessor<1, hp::DoFHandler<1> >::set_vertex_dof_index (const unsigned int /*vertex*/, - const unsigned int /*i*/, - const unsigned int /*index*/, - const unsigned int fe_index) const +void DoFObjectAccessor<1, hp::DoFHandler<1> >::set_vertex_dof_index (const unsigned int vertex, + const unsigned int i, + const unsigned int index, + const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<1>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - Assert (false, ExcInternalError()); + internal::hp::DoFLevel<0>::set_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i, + index); } template <> @@ -488,25 +490,29 @@ DoFObjectAccessor<1,hp::DoFHandler<1> >::set_dof_values (const Vector &/ #if deal_II_dimension == 2 template <> -void DoFObjectAccessor<1, hp::DoFHandler<2> >::set_vertex_dof_index (const unsigned int /*vertex*/, - const unsigned int /*i*/, - const unsigned int /*index*/, - const unsigned int fe_index) const +void DoFObjectAccessor<1, hp::DoFHandler<2> >::set_vertex_dof_index (const unsigned int vertex, + const unsigned int i, + const unsigned int index, + const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<2>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - Assert (false, ExcInternalError()); + internal::hp::DoFLevel<0>::set_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i, + index); } template <> -void DoFObjectAccessor<2, hp::DoFHandler<2> >::set_vertex_dof_index (const unsigned int /*vertex*/, - const unsigned int /*i*/, - const unsigned int /*index*/, - const unsigned int fe_index) const +void DoFObjectAccessor<2, hp::DoFHandler<2> >::set_vertex_dof_index (const unsigned int vertex, + const unsigned int i, + const unsigned int index, + const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<2>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - Assert (false, ExcInternalError()); + internal::hp::DoFLevel<0>::set_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i, + index); } template <> @@ -555,37 +561,43 @@ DoFObjectAccessor<2,hp::DoFHandler<2> >::set_dof_values (const Vector &/ #if deal_II_dimension == 3 template <> -void DoFObjectAccessor<1, hp::DoFHandler<3> >::set_vertex_dof_index (const unsigned int /*vertex*/, - const unsigned int /*i*/, - const unsigned int /*index*/, - const unsigned int fe_index) const +void DoFObjectAccessor<1, hp::DoFHandler<3> >::set_vertex_dof_index (const unsigned int vertex, + const unsigned int i, + const unsigned int index, + const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - Assert (false, ExcInternalError()); + internal::hp::DoFLevel<0>::set_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i, + index); } template <> -void DoFObjectAccessor<2, hp::DoFHandler<3> >::set_vertex_dof_index (const unsigned int /*vertex*/, - const unsigned int /*i*/, - const unsigned int /*index*/, - const unsigned int fe_index) const +void DoFObjectAccessor<2, hp::DoFHandler<3> >::set_vertex_dof_index (const unsigned int vertex, + const unsigned int i, + const unsigned int index, + const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - Assert (false, ExcInternalError()); + internal::hp::DoFLevel<0>::set_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i, + index); } template <> -void DoFObjectAccessor<3, hp::DoFHandler<3> >::set_vertex_dof_index (const unsigned int /*vertex*/, - const unsigned int /*i*/, - const unsigned int /*index*/, - const unsigned int fe_index) const +void DoFObjectAccessor<3, hp::DoFHandler<3> >::set_vertex_dof_index (const unsigned int vertex, + const unsigned int i, + const unsigned int index, + const unsigned int fe_index) const { - Assert (fe_index != hp::DoFHandler<3>::default_fe_index, - ExcMessage ("You need to specify a FE index when working with hp DoFHandlers")); - Assert (false, ExcInternalError()); + internal::hp::DoFLevel<0>::set_vertex_dof_index (*this->dof_handler, + this->vertex_index(vertex), + fe_index, + i, + index); } diff --git a/deal.II/deal.II/source/dofs/hp_dof_levels.cc b/deal.II/deal.II/source/dofs/hp_dof_levels.cc index c47a5dd720..a03584683b 100644 --- a/deal.II/deal.II/source/dofs/hp_dof_levels.cc +++ b/deal.II/deal.II/source/dofs/hp_dof_levels.cc @@ -55,10 +55,10 @@ namespace internal template unsigned int DoFLevel<0>:: - get_hp_vertex_dof_index (const ::hp::DoFHandler &dof_handler, - const unsigned int vertex_index, - const unsigned int fe_index, - const unsigned int local_index) + get_vertex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int vertex_index, + const unsigned int fe_index, + const unsigned int local_index) { Assert (fe_index != ::hp::DoFHandler::default_fe_index, ExcMessage ("You need to specify a FE index when working " @@ -99,11 +99,11 @@ namespace internal template void DoFLevel<0>:: - set_hp_vertex_dof_index (::hp::DoFHandler &dof_handler, - const unsigned int vertex_index, - const unsigned int fe_index, - const unsigned int local_index, - const unsigned int global_index) + set_vertex_dof_index (::hp::DoFHandler &dof_handler, + const unsigned int vertex_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index) { Assert (fe_index != ::hp::DoFHandler::default_fe_index, ExcMessage ("You need to specify a FE index when working " @@ -142,24 +142,477 @@ namespace internal } } + + + template + unsigned int + DoFLevel<1>:: + get_line_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int line_index, + const unsigned int fe_index, + const unsigned int local_index) const + { + Assert (fe_index != ::hp::DoFHandler::default_fe_index, + ExcMessage ("You need to specify a FE index when working " + "with hp DoFHandlers")); + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (local_index < dof_handler.get_fe()[fe_index].dofs_per_line, + ExcIndexRange(local_index, 0, + dof_handler.get_fe()[fe_index].dofs_per_line)); + + // if we are in 1d, then the + // only set of indices we store + // is the one for the cell, + // which is unique. then + // fe_index must be + // active_fe_index + if (dim == 1) + { + Assert (fe_index == this->active_fe_indices[line_index], + ExcMessage ("FE index does not match that of the present cell")); + return line_dofs[dof_line_index_offset[line_index]+local_index]; + } + else + { + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_line_index_offset[line_index]; + const unsigned int *pointer = &line_dofs[starting_offset]; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + if (*pointer == fe_index) + return *(pointer + 1 + local_index); + else + pointer += dof_handler.get_fe()[*pointer].dofs_per_line; + } + } + } + + + + template + void + DoFLevel<1>:: + set_line_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int line_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index) + { + Assert (fe_index != ::hp::DoFHandler::default_fe_index, + ExcMessage ("You need to specify a FE index when working " + "with hp DoFHandlers")); + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (local_index < dof_handler.get_fe()[fe_index].dofs_per_line, + ExcIndexRange(local_index, 0, + dof_handler.get_fe()[fe_index].dofs_per_line)); + + // if we are in 1d, then the + // only set of indices we store + // is the one for the cell, + // which is unique. then + // fe_index must be + // active_fe_index + if (dim == 1) + { + Assert (fe_index == this->active_fe_indices[line_index], + ExcMessage ("FE index does not match that of the present cell")); + line_dofs[dof_line_index_offset[line_index]+local_index] = global_index; + } + else + { + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_line_index_offset[line_index]; + unsigned int *pointer = &line_dofs[starting_offset]; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + if (*pointer == fe_index) + { + *(pointer + 1 + local_index) = global_index; + return; + } + else + pointer += dof_handler.get_fe()[*pointer].dofs_per_line; + } + } + } + + + + template + unsigned int + DoFLevel<2>:: + get_quad_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int quad_index, + const unsigned int fe_index, + const unsigned int local_index) const + { + Assert (dim >= 2, ExcMessage ("You can only access quads in 2d or higher")); + Assert (fe_index != ::hp::DoFHandler::default_fe_index, + ExcMessage ("You need to specify a FE index when working " + "with hp DoFHandlers")); + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (local_index < dof_handler.get_fe()[fe_index].dofs_per_quad, + ExcIndexRange(local_index, 0, + dof_handler.get_fe()[fe_index].dofs_per_quad)); + + // if we are in 1d, then the + // only set of indices we store + // is the one for the cell, + // which is unique. then + // fe_index must be + // active_fe_index + if (dim == 2) + { + Assert (fe_index == this->active_fe_indices[quad_index], + ExcMessage ("FE index does not match that of the present cell")); + return quad_dofs[dof_quad_index_offset[quad_index]+local_index]; + } + else + { + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_quad_index_offset[quad_index]; + const unsigned int *pointer = &quad_dofs[starting_offset]; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + if (*pointer == fe_index) + return *(pointer + 1 + local_index); + else + pointer += dof_handler.get_fe()[*pointer].dofs_per_quad; + } + } + } + + + + template + void + DoFLevel<2>:: + set_quad_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int quad_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index) + { + Assert (dim >= 2, ExcMessage ("You can only access quads in 2d or higher")); + Assert (fe_index != ::hp::DoFHandler::default_fe_index, + ExcMessage ("You need to specify a FE index when working " + "with hp DoFHandlers")); + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (local_index < dof_handler.get_fe()[fe_index].dofs_per_quad, + ExcIndexRange(local_index, 0, + dof_handler.get_fe()[fe_index].dofs_per_quad)); + + // if we are in 1d, then the + // only set of indices we store + // is the one for the cell, + // which is unique. then + // fe_index must be + // active_fe_index + if (dim == 2) + { + Assert (fe_index == this->active_fe_indices[quad_index], + ExcMessage ("FE index does not match that of the present cell")); + quad_dofs[dof_quad_index_offset[quad_index]+local_index] = global_index; + } + else + { + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_quad_index_offset[quad_index]; + unsigned int *pointer = &quad_dofs[starting_offset]; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + if (*pointer == fe_index) + { + *(pointer + 1 + local_index) = global_index; + return; + } + else + pointer += dof_handler.get_fe()[*pointer].dofs_per_quad; + } + } + } + + + + template + unsigned int + DoFLevel<3>:: + get_hex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int hex_index, + const unsigned int fe_index, + const unsigned int local_index) const + { + Assert (dim >= 3, ExcMessage ("You can only access hexs in 3d or higher")); + Assert (fe_index != ::hp::DoFHandler::default_fe_index, + ExcMessage ("You need to specify a FE index when working " + "with hp DoFHandlers")); + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (local_index < dof_handler.get_fe()[fe_index].dofs_per_hex, + ExcIndexRange(local_index, 0, + dof_handler.get_fe()[fe_index].dofs_per_hex)); + + // if we are in 1d, then the + // only set of indices we store + // is the one for the cell, + // which is unique. then + // fe_index must be + // active_fe_index + if (dim == 3) + { + Assert (fe_index == this->active_fe_indices[hex_index], + ExcMessage ("FE index does not match that of the present cell")); + return hex_dofs[dof_hex_index_offset[hex_index]+local_index]; + } + else + { + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_hex_index_offset[hex_index]; + const unsigned int *pointer = &hex_dofs[starting_offset]; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + if (*pointer == fe_index) + return *(pointer + 1 + local_index); + else + pointer += dof_handler.get_fe()[*pointer].dofs_per_hex; + } + } + } + + + + template + void + DoFLevel<3>:: + set_hex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int hex_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index) + { + Assert (dim >= 3, ExcMessage ("You can only access hexs in 3d or higher")); + Assert (fe_index != ::hp::DoFHandler::default_fe_index, + ExcMessage ("You need to specify a FE index when working " + "with hp DoFHandlers")); + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (local_index < dof_handler.get_fe()[fe_index].dofs_per_hex, + ExcIndexRange(local_index, 0, + dof_handler.get_fe()[fe_index].dofs_per_hex)); + + // if we are in 1d, then the + // only set of indices we store + // is the one for the cell, + // which is unique. then + // fe_index must be + // active_fe_index + if (dim == 3) + { + Assert (fe_index == this->active_fe_indices[hex_index], + ExcMessage ("FE index does not match that of the present cell")); + hex_dofs[dof_hex_index_offset[hex_index]+local_index] = global_index; + } + else + { + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_hex_index_offset[hex_index]; + unsigned int *pointer = &hex_dofs[starting_offset]; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + if (*pointer == fe_index) + { + *(pointer + 1 + local_index) = global_index; + return; + } + else + pointer += dof_handler.get_fe()[*pointer].dofs_per_hex; + } + } + } + // explicit instantiations template unsigned int DoFLevel<0>:: - get_hp_vertex_dof_index (const ::hp::DoFHandler &dof_handler, - const unsigned int vertex_index, - const unsigned int fe_index, - const unsigned int local_index); + get_vertex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int vertex_index, + const unsigned int fe_index, + const unsigned int local_index); template void DoFLevel<0>:: - set_hp_vertex_dof_index (::hp::DoFHandler &dof_handler, - const unsigned int vertex_index, - const unsigned int fe_index, - const unsigned int local_index, - const unsigned int global_index); + set_vertex_dof_index (::hp::DoFHandler &dof_handler, + const unsigned int vertex_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index); + + template + unsigned int + DoFLevel<1>:: + get_line_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int line_index, + const unsigned int fe_index, + const unsigned int local_index) const; + + template + void + DoFLevel<1>:: + set_line_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int line_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index); + +#if deal_II_dimension >= 2 + + template + unsigned int + DoFLevel<2>:: + get_quad_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int quad_index, + const unsigned int fe_index, + const unsigned int local_index) const; + + template + void + DoFLevel<2>:: + set_quad_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int quad_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index); + +#endif + +#if deal_II_dimension >= 3 + + template + unsigned int + DoFLevel<3>:: + get_hex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int hex_index, + const unsigned int fe_index, + const unsigned int local_index) const; + + template + void + DoFLevel<3>:: + set_hex_dof_index (const ::hp::DoFHandler &dof_handler, + const unsigned int hex_index, + const unsigned int fe_index, + const unsigned int local_index, + const unsigned int global_index); + +#endif } } -- 2.39.5