From 4026d0eca8a634be6ab426603d15a80386fc2b3a Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 2 May 2006 20:46:57 +0000 Subject: [PATCH] Implement reserve_space for all space dimensions. git-svn-id: https://svn.dealii.org/trunk@12975 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/hp_dof_handler.h | 9 + deal.II/deal.II/source/dofs/hp_dof_handler.cc | 712 ++++++++++++++---- 2 files changed, 560 insertions(+), 161 deletions(-) diff --git a/deal.II/deal.II/include/dofs/hp_dof_handler.h b/deal.II/deal.II/include/dofs/hp_dof_handler.h index a9655f818d..74b1c5d533 100644 --- a/deal.II/deal.II/include/dofs/hp_dof_handler.h +++ b/deal.II/deal.II/include/dofs/hp_dof_handler.h @@ -953,6 +953,15 @@ namespace hp */ void reserve_space (); + /** + * Do that part of reserving + * space that pertains to + * vertices, since this is the + * same in all space + * dimensions. + */ + void reserve_space_vertices (); + /** * Free all used memory. */ diff --git a/deal.II/deal.II/source/dofs/hp_dof_handler.cc b/deal.II/deal.II/source/dofs/hp_dof_handler.cc index 926ffa46e6..a33cb6a3c6 100644 --- a/deal.II/deal.II/source/dofs/hp_dof_handler.cc +++ b/deal.II/deal.II/source/dofs/hp_dof_handler.cc @@ -1576,7 +1576,7 @@ namespace hp DoFHandler<1>::distribute_dofs_on_cell (active_cell_iterator &cell, unsigned int next_free_dof) { - const FiniteElement<2> &fe = cell->get_fe(); + const FiniteElement<1> &fe = cell->get_fe(); const unsigned int fe_index = cell->active_fe_index (); // number dofs on vertices. to do @@ -1680,7 +1680,7 @@ namespace hp DoFHandler<3>::distribute_dofs_on_cell (active_cell_iterator &cell, unsigned int next_free_dof) { - const FiniteElement<2> &fe = cell->get_fe(); + const FiniteElement<3> &fe = cell->get_fe(); const unsigned int fe_index = cell->active_fe_index (); // number dofs on vertices. to do @@ -2058,16 +2058,13 @@ namespace hp template <> void DoFHandler<1>::reserve_space () { -//TODO[WB]: do the same as for 2d already + const unsigned int dim = 1; + Assert (finite_elements != 0, ExcNoFESelected()); Assert (finite_elements->size() > 0, ExcNoFESelected()); Assert (tria->n_levels() > 0, ExcInvalidTriangulation()); Assert (tria->n_levels() == levels.size (), ExcInternalError ()); - Assert (finite_elements->max_dofs_per_vertex() == 0, - ExcMessage ("hp finite elements are presently only supported " - "for discontinuous elements")); - // Release all space except the // active_fe_indices field which // we have to backup before @@ -2084,39 +2081,47 @@ namespace hp for (unsigned int level=0; leveln_levels(); ++level) { - levels.push_back (new internal::hp::DoFLevel<1>); - std::swap (active_fe_backup[level], levels[level]->active_fe_indices); + levels.push_back (new internal::hp::DoFLevel); + std::swap (active_fe_backup[level], + levels[level]->active_fe_indices); } } + // LINE (CELL) DOFs + // count how much space we need - // on each level and set the + // on each level for the cell + // dofs and set the // dof_*_offsets // data. initially set the latter // to an invalid index, and only // later set it to something // reasonable for active cells + // + // note that for cells, the + // situation is simpler than for + // other (lower dimensional) + // objects since exactly one + // finite element is used for it for (unsigned int level=0; leveln_levels(); ++level) { levels[level]->line_dof_offsets = std::vector (tria->n_raw_lines(level), invalid_dof_index); - unsigned int next_free_line_dof = 0; + unsigned int next_free_dof = 0; for (active_cell_iterator cell=begin_active(level); cell!=end_active(level); ++cell) if (!cell->has_children()) - { - levels[level]->line_dof_offsets[cell->index()] = next_free_line_dof; - next_free_line_dof += - (*finite_elements)[cell->active_fe_index()].dofs_per_line; - } + { + levels[level]->line_dof_offsets[cell->index()] = next_free_dof; + next_free_dof += cell->get_fe().dofs_per_line; + } - levels[level]->line_dofs = std::vector (next_free_line_dof, - invalid_dof_index); + levels[level]->line_dofs + = std::vector (next_free_dof, invalid_dof_index); } - // safety check: make sure that // the number of DoFs we // allocated is actually correct @@ -2127,14 +2132,13 @@ namespace hp #ifdef DEBUG for (unsigned int level=0; leveln_levels(); ++level) { - unsigned int n_line_dofs = 0; + unsigned int counter = 0; for (cell_iterator cell=begin_active(level); - cell!=end_active(level); ++cell) + cell!=end_active(level); ++cell) if (!cell->has_children()) - n_line_dofs += - (*finite_elements)[cell->active_fe_index()].dofs_per_line; + counter += cell->get_fe().dofs_per_line; - Assert (levels[level]->line_dofs.size() == n_line_dofs, ExcInternalError()); + Assert (levels[level]->line_dofs.size() == counter, ExcInternalError()); Assert (static_cast (std::count (levels[level]->line_dof_offsets.begin(), levels[level]->line_dof_offsets.end(), @@ -2144,9 +2148,12 @@ namespace hp ExcInternalError()); } #endif + + + // VERTEX DOFS + reserve_space_vertices (); } - #endif @@ -2178,42 +2185,46 @@ namespace hp for (unsigned int level=0; leveln_levels(); ++level) { - levels.push_back (new internal::hp::DoFLevel<2>); + levels.push_back (new internal::hp::DoFLevel); std::swap (active_fe_backup[level], levels[level]->active_fe_indices); } } + // QUAD (CELL) DOFs - // + // count how much space we need - // on each level for the quad + // on each level for the cell // dofs and set the // dof_*_offsets // data. initially set the latter // to an invalid index, and only // later set it to something // reasonable for active cells + // + // note that for cells, the + // situation is simpler than for + // other (lower dimensional) + // objects since exactly one + // finite element is used for it for (unsigned int level=0; leveln_levels(); ++level) { levels[level]->quad_dof_offsets = std::vector (tria->n_raw_quads(level), invalid_dof_index); - unsigned int next_free_quad_dof = 0; + unsigned int next_free_dof = 0; for (active_cell_iterator cell=begin_active(level); cell!=end_active(level); ++cell) if (!cell->has_children()) - { - levels[level]->quad_dof_offsets[cell->index()] - = next_free_quad_dof; - next_free_quad_dof - += (*finite_elements)[cell->active_fe_index()].dofs_per_quad; - } + { + levels[level]->quad_dof_offsets[cell->index()] = next_free_dof; + next_free_dof += cell->get_fe().dofs_per_quad; + } levels[level]->quad_dofs - = std::vector (next_free_quad_dof, - invalid_dof_index); + = std::vector (next_free_dof, invalid_dof_index); } // safety check: make sure that @@ -2226,15 +2237,13 @@ namespace hp #ifdef DEBUG for (unsigned int level=0; leveln_levels(); ++level) { - unsigned int n_quad_dofs = 0; + unsigned int counter = 0; for (cell_iterator cell=begin_active(level); - cell!=end_active(level); ++cell) + cell!=end_active(level); ++cell) if (!cell->has_children()) - n_quad_dofs += - (*finite_elements)[cell->active_fe_index()].dofs_per_quad; + counter += cell->get_fe().dofs_per_quad; - Assert (levels[level]->quad_dofs.size() == n_quad_dofs, - ExcInternalError()); + Assert (levels[level]->quad_dofs.size() == counter, ExcInternalError()); Assert (static_cast (std::count (levels[level]->quad_dof_offsets.begin(), levels[level]->quad_dof_offsets.end(), @@ -2244,7 +2253,7 @@ namespace hp ExcInternalError()); } #endif - + // LINE DOFS // @@ -2502,101 +2511,7 @@ namespace hp // VERTEX DOFS - - // The final step is to allocate - // vertex dof information. since - // vertices are sequentially - // numbered, what we do first is - // to set up an array in which we - // record whether a vertex is - // associated with any of the - // given fe's, by setting a - // bit. in a later step, we then - // actually allocate memory for - // the required dofs - { - std::vector > - vertex_fe_association (finite_elements->size(), - std::vector (tria->n_vertices(), false)); - - for (active_cell_iterator cell=begin_active(); cell!=end(); ++cell) - for (unsigned int v=0; v::vertices_per_cell; ++v) - vertex_fe_association[cell->active_fe_index()][cell->vertex_index(v)] - = true; - - // in debug mode, make sure - // that each vertex is - // associated with at least one - // fe (note that except for - // unused vertices, all - // vertices are actually - // active) -#ifdef DEBUG - for (unsigned int v=0; vn_vertices(); ++v) - if (tria->vertex_used(v) == true) - { - unsigned int fe=0; - for (; fesize(); ++fe) - if (vertex_fe_association[fe][v] == true) - break; - Assert (fe != finite_elements->size(), ExcInternalError()); - } -#endif - - // next count how much memory - // we actually need. for each - // vertex, we need one slot per - // fe to store the fe_index, - // plus dofs_per_vertex for - // this fe. in addition, we - // need one slot as the end - // marker for the - // fe_indices. at the same time - // already fill the - // vertex_dofs_offsets field - vertex_dofs_offsets.resize (tria->n_vertices(), - deal_II_numbers::invalid_unsigned_int); - - unsigned int vertex_slots_needed = 0; - for (unsigned int v=0; vn_vertices(); ++v) - if (tria->vertex_used(v) == true) - { - vertex_dofs_offsets[v] = vertex_slots_needed; - - for (unsigned int fe=0; fesize(); ++fe) - if (vertex_fe_association[fe][v] == true) - vertex_slots_needed += (*finite_elements)[fe].dofs_per_vertex + 1; - ++vertex_slots_needed; - } - - // now allocate the space we - // have determined we need, and - // set up the linked lists for - // each of the vertices - vertex_dofs.resize (vertex_slots_needed, invalid_dof_index); - for (unsigned int v=0; vn_vertices(); ++v) - if (tria->vertex_used(v) == true) - { - unsigned int pointer = vertex_dofs_offsets[v]; - for (unsigned int fe=0; fesize(); ++fe) - if (vertex_fe_association[fe][v] == true) - { - // if this vertex - // uses this fe, - // then set the - // fe_index and - // move the pointer - // ahead - vertex_dofs[pointer] = fe; - pointer += (*finite_elements)[fe].dofs_per_vertex + 1; - } - // finally place the end - // marker - vertex_dofs[pointer] = deal_II_numbers::invalid_unsigned_int; - ++pointer; - } - } - + reserve_space_vertices (); } #endif @@ -2607,16 +2522,13 @@ namespace hp template <> void DoFHandler<3>::reserve_space () { -//TODO[WB]: do the same as for 2d already + const unsigned int dim = 3; + Assert (finite_elements != 0, ExcNoFESelected()); Assert (finite_elements->size() > 0, ExcNoFESelected()); Assert (tria->n_levels() > 0, ExcInvalidTriangulation()); Assert (tria->n_levels() == levels.size (), ExcInternalError ()); - Assert (finite_elements->max_dofs_per_face() == 0, - ExcMessage ("hp finite elements are presently only supported " - "for discontinuous elements")); - // Release all space except the // active_fe_indices field which // we have to backup before @@ -2633,39 +2545,48 @@ namespace hp for (unsigned int level=0; leveln_levels(); ++level) { - levels.push_back (new internal::hp::DoFLevel<3>); - std::swap (active_fe_backup[level], levels[level]->active_fe_indices); + levels.push_back (new internal::hp::DoFLevel); + std::swap (active_fe_backup[level], + levels[level]->active_fe_indices); } } + + // HEX (CELL) DOFs + // count how much space we need - // on each level and set the + // on each level for the cell + // dofs and set the // dof_*_offsets // data. initially set the latter // to an invalid index, and only // later set it to something // reasonable for active cells + // + // note that for cells, the + // situation is simpler than for + // other (lower dimensional) + // objects since exactly one + // finite element is used for it for (unsigned int level=0; leveln_levels(); ++level) { levels[level]->hex_dof_offsets = std::vector (tria->n_raw_hexs(level), invalid_dof_index); - unsigned int next_free_hex_dof = 0; + unsigned int next_free_dof = 0; for (active_cell_iterator cell=begin_active(level); cell!=end_active(level); ++cell) if (!cell->has_children()) - { - levels[level]->hex_dof_offsets[cell->index()] = next_free_hex_dof; - next_free_hex_dof += - (*finite_elements)[cell->active_fe_index()].dofs_per_hex; - } + { + levels[level]->hex_dof_offsets[cell->index()] = next_free_dof; + next_free_dof += cell->get_fe().dofs_per_hex; + } - levels[level]->hex_dofs = std::vector (next_free_hex_dof, - invalid_dof_index); + levels[level]->hex_dofs + = std::vector (next_free_dof, invalid_dof_index); } - // safety check: make sure that // the number of DoFs we // allocated is actually correct @@ -2676,14 +2597,13 @@ namespace hp #ifdef DEBUG for (unsigned int level=0; leveln_levels(); ++level) { - unsigned int n_hex_dofs = 0; + unsigned int counter = 0; for (cell_iterator cell=begin_active(level); - cell!=end_active(level); ++cell) + cell!=end_active(level); ++cell) if (!cell->has_children()) - n_hex_dofs += - (*finite_elements)[cell->active_fe_index()].dofs_per_hex; + counter += cell->get_fe().dofs_per_hex; - Assert (levels[level]->hex_dofs.size() == n_hex_dofs, ExcInternalError()); + Assert (levels[level]->hex_dofs.size() == counter, ExcInternalError()); Assert (static_cast (std::count (levels[level]->hex_dof_offsets.begin(), levels[level]->hex_dof_offsets.end(), @@ -2693,12 +2613,482 @@ namespace hp ExcInternalError()); } #endif + + + // QUAD DOFS + // + // same here: count quad dofs, + // then allocate as much space as + // we need and prime the linked + // list for quad (see the + // description in hp::DoFLevels) + // with the indices we will + // need. note that our task is + // more complicated since two + // adjacent cells may have + // different active_fe_indices, + // in which case we need to + // allocate *two* sets of line + // dofs for the same line + // + // the way we do things is that + // we loop over all active cells + // (these are the ones that have + // DoFs only anyway) and all + // their faces. We note in the + // user flags whether we have + // previously visited a face and + // if so skip it (consequently, + // we have to save and later + // restore the line flags) + { + std::vector saved_quad_user_flags; + const_cast&>(*tria) + .save_user_flags_quad (saved_quad_user_flags); + const_cast&>(*tria).clear_user_flags_quad (); + + // an array to hold how many + // slots (see the hp::DoFLevel + // class) we will have to store + // on each level + std::vector n_quad_slots (tria->n_levels(), 0); + + for (active_cell_iterator cell=begin_active(); cell!=end(); ++cell) + for (unsigned int face=0; face::faces_per_cell; ++face) + if (! cell->face(face)->user_flag_set()) + { + // ok, face has not been + // visited. so we need to + // allocate space for + // it. let's see how much + // we need: we need one + // set if a) there is no + // neighbor behind this + // face, or b) the + // neighbor is not on the + // same level or further + // refined, or c) the + // neighbor is on the + // same level, but + // happens to have the + // same active_fe_index: + if (cell->at_boundary(face) + || + (cell->neighbor(face)->level() < cell->level()) + || + cell->neighbor(face)->has_children() + || + ((cell->neighbor(face)->level() == cell->level()) + && + !cell->neighbor(face)->has_children() + && + (cell->active_fe_index() == cell->neighbor(face)->active_fe_index()))) + // ok, one set of + // dofs. that makes + // one index, 1 times + // dofs_per_quad + // dofs, and one stop + // index + n_quad_slots[cell->level()] + += (*finite_elements)[cell->active_fe_index()].dofs_per_quad + 2; + + // otherwise we do + // indeed need two + // sets, i.e. two + // indices, two sets of + // dofs, and one stop + // index: + else + n_quad_slots[cell->level()] + += ((*finite_elements)[cell->active_fe_index()].dofs_per_quad + + + (*finite_elements)[cell->neighbor(face)->active_fe_index()] + .dofs_per_quad + + + 3); + + // mark this face as + // visited + cell->face(face)->set_user_flag (); + } + + // now that we know how many + // quad dofs we will have to + // have on each level, allocate + // the memory. note that we + // allocate offsets for all + // quads, though only the + // active ones will have a + // non-invalid value later on + for (unsigned int level=0; leveln_levels(); ++level) + { + levels[level]->quad_dof_offsets + = std::vector (tria->n_raw_quads(level), + invalid_dof_index); + levels[level]->quad_dofs + = std::vector (n_quad_slots[level], + invalid_dof_index); + } + + // with the memory now + // allocated, loop over the + // cells again and prime the + // _offset values as well as + // the fe_index fields + const_cast&>(*tria).clear_user_flags_quad (); + + std::vector next_free_quad_slot (tria->n_levels(), 0); + + for (active_cell_iterator cell=begin_active(); cell!=end(); ++cell) + for (unsigned int face=0; face::faces_per_cell; ++face) + if (! cell->face(face)->user_flag_set()) + { + // same decision tree + // as before + if (cell->at_boundary(face) + || + (cell->neighbor(face)->level() < cell->level()) + || + cell->neighbor(face)->has_children() + || + ((cell->neighbor(face)->level() == cell->level()) + && + !cell->neighbor(face)->has_children() + && + (cell->active_fe_index() == cell->neighbor(face)->active_fe_index()))) + { + levels[cell->level()] + ->quad_dof_offsets[cell->face(face)->index()] + = next_free_quad_slot[cell->level()]; + + // set first slot + // for this quad to + // active_fe_index + // of this face + levels[cell->level()] + ->quad_dofs[next_free_quad_slot[cell->level()]] + = cell->active_fe_index(); + + // the next + // dofs_per_quad + // indices remain + // unset for the + // moment (i.e. at + // invalid_dof_index). + // following this + // comes the stop + // index, which + // also is + // invalid_dof_index + // and therefore + // does not have to + // be explicitly + // set + + // finally, mark + // those slots as + // used + next_free_quad_slot[cell->level()] + += (*finite_elements)[cell->active_fe_index()].dofs_per_quad + 2; + } + else + { + levels[cell->level()] + ->quad_dof_offsets[cell->face(face)->index()] + = next_free_quad_slot[cell->level()]; + + // set first slot + // for this quad to + // active_fe_index + // of this face + levels[cell->level()] + ->quad_dofs[next_free_quad_slot[cell->level()]] + = cell->active_fe_index(); + + // the next + // dofs_per_quad + // indices remain + // unset for the + // moment (i.e. at + // invalid_dof_index). + // + // then comes the + // fe_index for the + // neighboring + // cell: + levels[cell->level()] + ->quad_dofs[next_free_quad_slot[cell->level()] + + + (*finite_elements)[cell->active_fe_index()].dofs_per_quad + + + 1] + = cell->neighbor(face)->active_fe_index(); + // then again a set + // of dofs that we + // need not set + // right now + // + // following this + // comes the stop + // index, which + // also is + // invalid_dof_index + // and therefore + // does not have to + // be explicitly + // set + + // finally, mark + // those slots as + // used + next_free_quad_slot[cell->level()] + += ((*finite_elements)[cell->active_fe_index()].dofs_per_quad + + + (*finite_elements)[cell->neighbor(face)->active_fe_index()] + .dofs_per_quad + + + 3); + } + + // mark this face as + // visited + cell->face(face)->set_user_flag (); + } + + // we should have moved the + // cursor for each level to the + // total number of dofs on that + // level. check that + for (unsigned int level=0; leveln_levels(); ++level) + Assert (next_free_quad_slot[level] == n_quad_slots[level], + ExcInternalError()); + + // at the end, restore the user + // flags for the quads + const_cast&>(*tria) + .load_user_flags_quad (saved_quad_user_flags); + } + + + // LINE DOFS + + // the situation here is pretty + // much like with vertices: there + // can be an arbitrary number of + // finite elements associated + // with each line. the situation + // is more complicated, however, + // since lines have no global + // ordering, but are rather + // organized into levels. + // + // the algorithm we use is + // somewhat similar to what we do + // in reserve_space_vertices(), + // except that we have to to + // separate work for all levels + for (unsigned int level=0; leveln_levels(); ++level) + { + // what we do first is to set up + // an array in which we record + // whether a line is associated + // with any of the given fe's, by + // setting a bit. in a later + // step, we then actually + // allocate memory for the + // required dofs + std::vector > + line_fe_association (finite_elements->size(), + std::vector (tria->n_raw_lines(level), + false)); + + for (active_cell_iterator cell=begin_active(level); + cell!=end_active(level); ++cell) + for (unsigned int l=0; l::lines_per_cell; ++l) + line_fe_association[cell->active_fe_index()][cell->line_index(l)] + = true; + + // first check which of the + // lines is used at all, + // i.e. is associated with a + // finite element. we do this + // since not all lines may + // actually be used, in which + // case we do not have to + // allocate any memory at + // all + std::vector line_is_used (tria->n_raw_lines(level), false); + for (unsigned int line=0; linen_raw_lines(level); ++line) + for (unsigned int fe=0; fesize(); ++fe) + if (line_fe_association[fe][line] == true) + { + line_is_used[line] = true; + break; + } + + // next count how much memory + // we actually need. for each + // line, we need one slot per + // fe to store the fe_index, + // plus dofs_per_line for + // this fe. in addition, we + // need one slot as the end + // marker for the + // fe_indices. at the same + // time already fill the + // line_dofs_offsets field + levels[level]->line_dof_offsets + .resize (tria->n_raw_lines(level), + deal_II_numbers::invalid_unsigned_int); + + unsigned int line_slots_needed = 0; + for (unsigned int line=0; linen_raw_lines(level); ++line) + if (line_is_used[line] == true) + { + levels[level]->line_dof_offsets[line] = line_slots_needed; + + for (unsigned int fe=0; fesize(); ++fe) + if (line_fe_association[fe][line] == true) + line_slots_needed += (*finite_elements)[fe].dofs_per_line + 1; + ++line_slots_needed; + } + + // now allocate the space we + // have determined we need, and + // set up the linked lists for + // each of the lines + levels[level]->line_dofs.resize (line_slots_needed, invalid_dof_index); + for (unsigned int line=0; linen_raw_lines(level); ++line) + if (line_is_used[line] == true) + { + unsigned int pointer = levels[level]->line_dof_offsets[line]; + for (unsigned int fe=0; fesize(); ++fe) + if (line_fe_association[fe][line] == true) + { + // if this line + // uses this fe, + // then set the + // fe_index and + // move the + // pointer ahead + levels[level]->line_dofs[pointer] = fe; + pointer += (*finite_elements)[fe].dofs_per_line + 1; + } + // finally place the end + // marker + levels[level]->line_dofs[pointer] = deal_II_numbers::invalid_unsigned_int; + } + } + + + // VERTEX DOFS + reserve_space_vertices (); } #endif + + template + void + DoFHandler::reserve_space_vertices () + { + // The final step is allocating + // memory is to set up vertex dof + // information. since vertices + // are sequentially numbered, + // what we do first is to set up + // an array in which we record + // whether a vertex is associated + // with any of the given fe's, by + // setting a bit. in a later + // step, we then actually + // allocate memory for the + // required dofs + std::vector > + vertex_fe_association (finite_elements->size(), + std::vector (tria->n_vertices(), false)); + + for (active_cell_iterator cell=begin_active(); cell!=end(); ++cell) + for (unsigned int v=0; v::vertices_per_cell; ++v) + vertex_fe_association[cell->active_fe_index()][cell->vertex_index(v)] + = true; + + // in debug mode, make sure + // that each vertex is + // associated with at least one + // fe (note that except for + // unused vertices, all + // vertices are actually + // active) +#ifdef DEBUG + for (unsigned int v=0; vn_vertices(); ++v) + if (tria->vertex_used(v) == true) + { + unsigned int fe=0; + for (; fesize(); ++fe) + if (vertex_fe_association[fe][v] == true) + break; + Assert (fe != finite_elements->size(), ExcInternalError()); + } +#endif + + // next count how much memory + // we actually need. for each + // vertex, we need one slot per + // fe to store the fe_index, + // plus dofs_per_vertex for + // this fe. in addition, we + // need one slot as the end + // marker for the + // fe_indices. at the same time + // already fill the + // vertex_dofs_offsets field + vertex_dofs_offsets.resize (tria->n_vertices(), + deal_II_numbers::invalid_unsigned_int); + + unsigned int vertex_slots_needed = 0; + for (unsigned int v=0; vn_vertices(); ++v) + if (tria->vertex_used(v) == true) + { + vertex_dofs_offsets[v] = vertex_slots_needed; + + for (unsigned int fe=0; fesize(); ++fe) + if (vertex_fe_association[fe][v] == true) + vertex_slots_needed += (*finite_elements)[fe].dofs_per_vertex + 1; + ++vertex_slots_needed; + } + + // now allocate the space we + // have determined we need, and + // set up the linked lists for + // each of the vertices + vertex_dofs.resize (vertex_slots_needed, invalid_dof_index); + for (unsigned int v=0; vn_vertices(); ++v) + if (tria->vertex_used(v) == true) + { + unsigned int pointer = vertex_dofs_offsets[v]; + for (unsigned int fe=0; fesize(); ++fe) + if (vertex_fe_association[fe][v] == true) + { + // if this vertex + // uses this fe, + // then set the + // fe_index and + // move the pointer + // ahead + vertex_dofs[pointer] = fe; + pointer += (*finite_elements)[fe].dofs_per_vertex + 1; + } + // finally place the end + // marker + vertex_dofs[pointer] = deal_II_numbers::invalid_unsigned_int; + } + } + + + template void DoFHandler::create_active_fe_table () { -- 2.39.5